Coverage for biobb_vs/gnina/common.py: 71%
112 statements
« prev ^ index » next coverage.py v7.15.3, created at 2026-08-03 13:34 +0000
1"""Common functions for package biobb_vs.gnina"""
3import gzip
4import json
5import re
6from pathlib import Path, PurePath
8from biobb_common.tools import file_utils as fu
10# Matches an SDF data field header, e.g. "> <minimizedAffinity>"
11SDF_TAG_PATTERN = re.compile(r">\s*<([^>]+)>")
12# Record terminator of the SDF format
13SDF_TERMINATOR = "$$$$"
16# CHECK PARAMETERS
19def check_input_path(path, argument, out_log, classname):
20 """Checks input file"""
21 if not Path(path).exists():
22 fu.log(classname + ': Unexisting %s file, exiting' % argument, out_log)
23 raise SystemExit(classname + ': Unexisting %s file' % argument)
24 file_extension = PurePath(path).suffix
25 if not is_valid_file(file_extension[1:], argument):
26 fu.log(classname + ': Format %s in %s file is not compatible' % (file_extension[1:], argument), out_log)
27 raise SystemExit(classname + ': Format %s in %s file is not compatible' % (file_extension[1:], argument))
28 check_gzip_extension(path, argument, out_log, classname)
29 return path
32def check_output_path(path, argument, optional, out_log, classname):
33 """Checks output file"""
34 if optional and not path:
35 return None
36 if PurePath(path).parent and not Path(PurePath(path).parent).exists():
37 fu.log(classname + ': Unexisting %s folder, exiting' % argument, out_log)
38 raise SystemExit(classname + ': Unexisting %s folder' % argument)
39 file_extension = PurePath(path).suffix
40 if not is_valid_file(file_extension[1:], argument):
41 fu.log(classname + ': Format %s in %s file is not compatible' % (file_extension[1:], argument), out_log)
42 raise SystemExit(classname + ': Format %s in %s file is not compatible' % (file_extension[1:], argument))
43 check_gzip_extension(path, argument, out_log, classname)
44 return path
47def is_valid_file(ext, argument):
48 """Checks if file format is compatible"""
49 formats = {
50 'input_ligand_path': ['sdf', 'mol2', 'pdb', 'pdbqt'],
51 'input_receptor_path': ['pdb', 'pdbqt'],
52 'input_box_path': ['pdb'],
53 'input_autobox_path': ['sdf', 'mol2', 'pdb', 'pdbqt', 'pqr'],
54 'input_sdf_path': ['sdf', 'gz'],
55 'output_sdf_path': ['sdf', 'gz'],
56 'output_summary_path': ['json'],
57 'output_log_path': ['log']
58 }
59 return ext in formats[argument]
62def check_gzip_extension(path, argument, out_log, classname):
63 """Checks that a gzip compressed file also declares the format it wraps
65 gnina and Open Babel take the molecular format from the file extension, and
66 a bare .gz gives them nothing to work with. Only .sdf.gz is accepted.
67 """
68 suffixes = PurePath(path).suffixes
69 if suffixes[-1:] == ['.gz'] and suffixes[-2:-1] != ['.sdf']:
70 fu.log(classname + ': %s must use a .sdf.gz extension to be gzip compressed' % argument, out_log)
71 raise SystemExit(classname + ': %s must use a .sdf.gz extension to be gzip compressed' % argument)
72 return path
75# READ / WRITE SDF
78def open_sdf(path, mode='rt'):
79 """Opens a plain or gzip compressed SDF file, transparently"""
80 if PurePath(path).suffix == '.gz':
81 return gzip.open(path, mode, newline='')
82 return open(path, mode, newline='')
85def read_sdf_records(path):
86 """Splits a multi record SDF file into its records
88 Returns a list of dictionaries, one per record, each holding:
89 * **text**: the raw record text, terminator included, so that a record
90 can be written back out verbatim.
91 * **name**: the molecule name, that is the first line of the record.
92 * **data**: the SD data fields of the record, as a tag to value mapping.
93 * **ligand_index**: 1 based index of the ligand the record belongs to.
94 * **pose**: 1 based rank of the record within its ligand.
96 gnina writes the poses of every input ligand consecutively and does not tag
97 them with a ligand identifier, so records are grouped into ligands whenever
98 the molecule name changes. Consecutive ligands sharing the same name are
99 therefore reported as a single ligand.
100 """
101 with open_sdf(path, 'rt') as sdf_file:
102 content = sdf_file.read()
104 records = []
105 record_lines = []
106 for line in content.splitlines(keepends=True):
107 record_lines.append(line)
108 if line.strip() == SDF_TERMINATOR:
109 records.append(parse_sdf_record(record_lines))
110 record_lines = []
112 # keep a trailing record that is missing its terminator instead of losing it
113 if any(line.strip() for line in record_lines):
114 records.append(parse_sdf_record(record_lines))
116 # group the records into ligands, gnina writes each ligand's poses together
117 ligand_index = 0
118 previous_name = None
119 pose = 0
120 for record in records:
121 if record['name'] != previous_name:
122 ligand_index += 1
123 previous_name = record['name']
124 pose = 0
125 pose += 1
126 record['ligand_index'] = ligand_index
127 record['pose'] = pose
129 return records
132def parse_sdf_record(record_lines):
133 """Parses a single SDF record into its name, data fields and raw text"""
134 name = record_lines[0].strip() if record_lines else ''
135 data = {}
136 tag = None
137 values = []
139 for line in record_lines:
140 if line.strip() == SDF_TERMINATOR:
141 break
142 match = SDF_TAG_PATTERN.match(line)
143 if match:
144 if tag is not None:
145 data[tag] = '\n'.join(values)
146 tag = match.group(1)
147 values = []
148 continue
149 if tag is not None:
150 if not line.strip():
151 data[tag] = '\n'.join(values)
152 tag = None
153 values = []
154 else:
155 values.append(line.strip())
157 if tag is not None:
158 data[tag] = '\n'.join(values)
160 return {'text': ''.join(record_lines), 'name': name, 'data': data}
163def to_number(value):
164 """Converts an SD data value to a float, leaving it untouched if it is not numeric"""
165 try:
166 return float(value)
167 except (TypeError, ValueError):
168 return value
171# PROCESS OUTPUTS
174def process_output_gnina(output_sdf_path, output_summary_path, out_log, classname):
175 """Generates the output_summary_path JSON file from the poses written by gnina
177 Every pose becomes one entry holding the ligand it belongs to, its rank
178 within that ligand and every score gnina attached to it as an SD data field
179 (minimizedAffinity, CNNscore, CNNaffinity, CNN_VS and, for model ensembles,
180 CNNvariance). Which scores are present depends on the cnn_scoring property.
181 """
182 if not Path(output_sdf_path).exists():
183 fu.log(classname + ': %s not found, skipping the summary' % output_sdf_path, out_log)
184 raise SystemExit(classname + ': Error executing gnina, %s was not created' % output_sdf_path)
186 records = read_sdf_records(output_sdf_path)
188 data = {}
189 for index, record in enumerate(records, start=1):
190 entry = {
191 'ligand_name': record['name'],
192 'ligand_index': record['ligand_index'],
193 'pose': record['pose']
194 }
195 for tag, value in record['data'].items():
196 entry[tag] = to_number(value)
197 data['pose' + str(index)] = entry
199 ligands = len({record['ligand_index'] for record in records})
200 fu.log('%d poses found for %d ligand(s)' % (len(records), ligands), out_log)
202 fu.log('Saving summary to %s file' % output_summary_path, out_log)
203 with open(output_summary_path, 'w') as outfile:
204 json.dump(data, outfile, indent=4)
206 return data