Mercurial > repos > bgruening > openbabel_change_title
comparison distance_finder.py @ 16:a638d8d13bb3 draft default tip
planemo upload for repository https://github.com/bgruening/galaxytools/tree/master/chemicaltoolbox/openbabel commit d9c51279c061a1da948a2582d5b502ca7573adbf
| author | bgruening |
|---|---|
| date | Thu, 15 Aug 2024 11:02:57 +0000 |
| parents | 2cd8aee0d830 |
| children |
comparison
equal
deleted
inserted
replaced
| 15:01ab890bf0ad | 16:a638d8d13bb3 |
|---|---|
| 17 | 17 |
| 18 from openbabel import pybel | 18 from openbabel import pybel |
| 19 | 19 |
| 20 | 20 |
| 21 def log(*args, **kwargs): | 21 def log(*args, **kwargs): |
| 22 """Log output to STDERR | 22 """Log output to STDERR""" |
| 23 """ | 23 print(*args, file=sys.stderr, **kwargs) |
| 24 print(*args, file=sys.stderr, ** kwargs) | |
| 25 | 24 |
| 26 | 25 |
| 27 def execute(ligands_sdf, points_file, outfile): | 26 def execute(ligands_sdf, points_file, outfile): |
| 28 """ | 27 """ |
| 29 :param ligands_sdf: A SDF with the 3D molecules to test | 28 :param ligands_sdf: A SDF with the 3D molecules to test |
| 33 """ | 32 """ |
| 34 | 33 |
| 35 points = [] | 34 points = [] |
| 36 | 35 |
| 37 # read the points | 36 # read the points |
| 38 with open(points_file, 'r') as f: | 37 with open(points_file, "r") as f: |
| 39 for line in f.readlines(): | 38 for line in f.readlines(): |
| 40 line.strip() | 39 line.strip() |
| 41 if line: | 40 if line: |
| 42 p = line.split() | 41 p = line.split() |
| 43 if len(p) == 3: | 42 if len(p) == 3: |
| 44 points.append((float(p[0]), float(p[1]), float(p[2]))) | 43 points.append((float(p[0]), float(p[1]), float(p[2]))) |
| 45 log("Read points", p) | 44 log("Read points", p) |
| 46 continue | 45 continue |
| 47 log("Failed to read line:", line) | 46 log("Failed to read line:", line) |
| 48 log('Found', len(points), 'atom points') | 47 log("Found", len(points), "atom points") |
| 49 | 48 |
| 50 sdf_writer = pybel.Outputfile("sdf", outfile, overwrite=True) | 49 sdf_writer = pybel.Outputfile("sdf", outfile, overwrite=True) |
| 51 | 50 |
| 52 count = 0 | 51 count = 0 |
| 53 for mol in pybel.readfile("sdf", ligands_sdf): | 52 for mol in pybel.readfile("sdf", ligands_sdf): |
| 54 count += 1 | 53 count += 1 |
| 55 if count % 50000 == 0: | 54 if count % 50000 == 0: |
| 56 log('Processed', count) | 55 log("Processed", count) |
| 57 | 56 |
| 58 try: | 57 try: |
| 59 # print("Processing mol", mol.title) | 58 # print("Processing mol", mol.title) |
| 60 clone = pybel.Molecule(mol) | 59 clone = pybel.Molecule(mol) |
| 61 clone.removeh() | 60 clone.removeh() |
| 68 for point in points: | 67 for point in points: |
| 69 p += 1 | 68 p += 1 |
| 70 distances = [] | 69 distances = [] |
| 71 for i in coords: | 70 for i in coords: |
| 72 # calculates distance based on cartesian coordinates | 71 # calculates distance based on cartesian coordinates |
| 73 distance = math.sqrt((point[0] - i[0])**2 + (point[1] - i[1])**2 + (point[2] - i[2])**2) | 72 distance = math.sqrt( |
| 73 (point[0] - i[0]) ** 2 | |
| 74 + (point[1] - i[1]) ** 2 | |
| 75 + (point[2] - i[2]) ** 2 | |
| 76 ) | |
| 74 distances.append(distance) | 77 distances.append(distance) |
| 75 # log("distance:", distance) | 78 # log("distance:", distance) |
| 76 min_distance = min(distances) | 79 min_distance = min(distances) |
| 77 # log('Min:', min_distance) | 80 # log('Min:', min_distance) |
| 78 # log(count, p, min_distance) | 81 # log(count, p, min_distance) |
| 79 | 82 |
| 80 mol.data['distance' + str(p)] = min_distance | 83 mol.data["distance" + str(p)] = min_distance |
| 81 | 84 |
| 82 sdf_writer.write(mol) | 85 sdf_writer.write(mol) |
| 83 | 86 |
| 84 except Exception as e: | 87 except Exception as e: |
| 85 log('Failed to handle molecule: ' + str(e)) | 88 log("Failed to handle molecule: " + str(e)) |
| 86 continue | 89 continue |
| 87 | 90 |
| 88 sdf_writer.close() | 91 sdf_writer.close() |
| 89 log('Wrote', count, 'molecules') | 92 log("Wrote", count, "molecules") |
| 90 | 93 |
| 91 | 94 |
| 92 def main(): | 95 def main(): |
| 93 global work_dir | 96 global work_dir |
| 94 | 97 |
| 95 parser = argparse.ArgumentParser(description='XChem distances - measure distances to particular points') | 98 parser = argparse.ArgumentParser( |
| 96 parser.add_argument('-i', '--input', help="SDF containing the 3D molecules to score)") | 99 description="XChem distances - measure distances to particular points" |
| 97 parser.add_argument('-p', '--points', help="PDB format file with atoms") | 100 ) |
| 98 parser.add_argument('-o', '--outfile', default='output.sdf', help="File name for results") | 101 parser.add_argument( |
| 102 "-i", "--input", help="SDF containing the 3D molecules to score)" | |
| 103 ) | |
| 104 parser.add_argument("-p", "--points", help="PDB format file with atoms") | |
| 105 parser.add_argument( | |
| 106 "-o", "--outfile", default="output.sdf", help="File name for results" | |
| 107 ) | |
| 99 | 108 |
| 100 args = parser.parse_args() | 109 args = parser.parse_args() |
| 101 log("XChem distances args: ", args) | 110 log("XChem distances args: ", args) |
| 102 | 111 |
| 103 execute(args.input, args.points, args.outfile) | 112 execute(args.input, args.points, args.outfile) |
