Coverage for biobb_gromacs/gromacs_extra/ndx2resttop.py: 68%
120 statements
« prev ^ index » next coverage.py v7.9.1, created at 2025-06-25 09:23 +0000
« prev ^ index » next coverage.py v7.9.1, created at 2025-06-25 09:23 +0000
1#!/usr/bin/env python3
3"""Module containing the Ndx2resttop class and the command line interface."""
4import fnmatch
5import argparse
6from typing import Optional
7from typing import Any
8from pathlib import Path
9from biobb_common.generic.biobb_object import BiobbObject
10from biobb_common.configuration import settings
11from biobb_common.tools import file_utils as fu
12from biobb_common.tools.file_utils import launchlogger
15class Ndx2resttop(BiobbObject):
16 """
17 | biobb_gromacs Ndx2resttop
18 | Generate a restrained topology from an index NDX file.
19 | This module automatizes the process of restrained topology generation starting from an index NDX file.
21 Args:
22 input_ndx_path (str): Path to the input NDX index file. File type: input. `Sample file <https://github.com/bioexcel/biobb_gromacs/raw/master/biobb_gromacs/test/data/gromacs_extra/ndx2resttop.ndx>`_. Accepted formats: ndx (edam:format_2033).
23 input_top_zip_path (str): Path the input TOP topology in zip format. File type: input. `Sample file <https://github.com/bioexcel/biobb_gromacs/raw/master/biobb_gromacs/test/data/gromacs_extra/ndx2resttop.zip>`_. Accepted formats: zip (edam:format_3987).
24 output_top_zip_path (str): Path the output TOP topology in zip format. File type: output. `Sample file <https://github.com/bioexcel/biobb_gromacs/raw/master/biobb_gromacs/test/reference/gromacs_extra/ref_ndx2resttop.zip>`_. Accepted formats: zip (edam:format_3987).
25 properties (dict - Python dictionary object containing the tool parameters, not input/output files):
26 * **posres_names** (*list*) - ("CUSTOM_POSRES") String with names of the position restraints to be included in the topology file separated by spaces. If provided it should match the length of the ref_rest_chain_triplet_list.
27 * **force_constants** (*str*) - ("500 500 500") Array of three floats defining the force constants.
28 * **ref_rest_chain_triplet_list** (*str*) - (None) Triplet list composed by (reference group, restrain group, chain) list.
30 Examples:
31 This is a use example of how to use the building block from Python::
33 from biobb_gromacs.gromacs_extra.ndx2resttop import ndx2resttop
34 prop = { 'ref_rest_chain_triplet_list': '( Chain_A, Chain_A_noMut, A ), ( Chain_B, Chain_B_noMut, B ), ( Chain_C, Chain_C_noMut, C ), ( Chain_D, Chain_D_noMut, D )' }
35 ndx2resttop(input_ndx_path='/path/to/myIndex.ndx',
36 input_top_zip_path='/path/to/myTopology.zip',
37 output_top_zip_path='/path/to/newTopology.zip',
38 properties=prop)
40 Info:
41 * wrapped_software:
42 * name: In house
43 * license: Apache-2.0
44 * ontology:
45 * name: EDAM
46 * schema: http://edamontology.org/EDAM.owl
47 """
49 def __init__(self, input_ndx_path: str, input_top_zip_path: str, output_top_zip_path: str,
50 properties: Optional[dict] = None, **kwargs) -> None:
51 properties = properties or {}
53 # Call parent class constructor
54 super().__init__(properties)
55 self.locals_var_dict = locals().copy()
57 # Input/Output files
58 self.io_dict = {
59 "in": {"input_ndx_path": input_ndx_path, "input_top_zip_path": input_top_zip_path},
60 "out": {"output_top_zip_path": output_top_zip_path}
61 }
63 # Properties specific for BB
64 self.posres_names = properties.get('posres_names')
65 self.force_constants = properties.get('force_constants', '500 500 500')
66 self.ref_rest_chain_triplet_list = properties.get('ref_rest_chain_triplet_list')
68 # Check the properties
69 self.check_properties(properties)
70 self.check_arguments()
72 @launchlogger
73 def launch(self) -> int:
74 """Execute the :class:`Ndx2resttop <gromacs_extra.ndx2resttop.Ndx2resttop>` object."""
75 # Setup Biobb
76 if self.check_restart():
77 return 0
79 top_file = fu.unzip_top(zip_file=self.io_dict['in'].get("input_top_zip_path", ""), out_log=self.out_log)
81 # Create index list of index file (dictionary with the index of the start and stop lines of each group) :)
82 index_dic: dict[str, Any] = {}
83 lines = open(self.io_dict['in'].get("input_ndx_path", "")).read().splitlines()
84 for index, line in enumerate(lines):
86 # New group
87 if line.startswith('['):
88 index_dic[line] = [index, 0]
90 # Update current group
91 label = line
93 # Close previous group
94 if index > 0:
95 index_dic[label] = [index_dic[label][0], index]
97 # Last group of the file
98 if index == len(lines)-1:
99 index_dic[label] = [index_dic[label][0], index]
101 # Catch groups with just one line
102 for label in index_dic.keys():
103 if (index_dic[label][0]+1) == index_dic[label][1]:
104 index_dic[label][1] += 1
106 fu.log('Index_dic: '+str(index_dic), self.out_log, self.global_log)
108 self.ref_rest_chain_triplet_list = [tuple(elem.strip(' ()').replace(' ', '').split(',')) for elem in str(self.ref_rest_chain_triplet_list).split('),')]
109 fu.log('ref_rest_chain_triplet_list: ' + str(self.ref_rest_chain_triplet_list), self.out_log, self.global_log)
111 if self.posres_names:
112 self.posres_names = [elem.strip() for elem in self.posres_names.split()]
113 fu.log('posres_names: ' + str(self.posres_names), self.out_log, self.global_log)
114 else:
115 self.posres_names = ['CUSTOM_POSRES']*len(self.ref_rest_chain_triplet_list)
116 fu.log('posres_names: ' + str(self.posres_names), self.out_log, self.global_log)
118 # Check if the number of posres_names matches the number of ref_rest_chain_triplet_list
119 if len(self.posres_names) != len(self.ref_rest_chain_triplet_list):
120 raise ValueError("If posres_names is provided, it should match the number of ref_rest_chain_triplet_list")
122 for triplet, posre_name in zip(self.ref_rest_chain_triplet_list, self.posres_names):
124 reference_group, restrain_group, chain = triplet
126 fu.log('Reference group: '+reference_group, self.out_log, self.global_log)
127 fu.log('Restrain group: '+restrain_group, self.out_log, self.global_log)
128 fu.log('Chain: '+chain, self.out_log, self.global_log)
129 self.io_dict['out']["output_itp_path"] = fu.create_name(path=str(Path(top_file).parent), prefix=self.prefix, step=self.step, name=restrain_group+'_posre.itp')
131 # Mapping atoms from absolute enumeration to Chain relative enumeration
132 fu.log('reference_group_index: start_closed:'+str(index_dic['[ '+reference_group+' ]'][0]+1)+' stop_open: '+str(index_dic['[ '+reference_group+' ]'][1]), self.out_log, self.global_log)
133 reference_group_list = [int(elem) for line in lines[index_dic['[ '+reference_group+' ]'][0]+1: index_dic['[ '+reference_group+' ]'][1]] for elem in line.split()]
134 fu.log('restrain_group_index: start_closed:'+str(index_dic['[ '+restrain_group+' ]'][0]+1)+' stop_open: '+str(index_dic['[ '+restrain_group+' ]'][1]), self.out_log, self.global_log)
135 restrain_group_list = [int(elem) for line in lines[index_dic['[ '+restrain_group+' ]'][0]+1: index_dic['[ '+restrain_group+' ]'][1]] for elem in line.split()]
136 selected_list = [reference_group_list.index(atom)+1 for atom in restrain_group_list]
137 # Creating new ITP with restrains
138 with open(self.io_dict['out'].get("output_itp_path", ''), 'w') as f:
139 fu.log('Creating: '+str(f)+' and adding the selected atoms force constants', self.out_log, self.global_log)
140 f.write('[ position_restraints ]\n')
141 f.write('; atom type fx fy fz\n')
142 for atom in selected_list:
143 f.write(str(atom)+' 1 '+self.force_constants+'\n')
145 multi_chain = False
146 # For multi-chain topologies
147 # Including new ITP in the corresponding ITP-chain file
148 for file_dir in Path(top_file).parent.iterdir():
149 if "posre" not in file_dir.name and not file_dir.name.endswith("_pr.itp"):
150 if fnmatch.fnmatch(str(file_dir), "*_chain_"+chain+".itp"):
151 multi_chain = True
152 with open(str(file_dir), 'a') as f:
153 fu.log('Opening: '+str(f)+' and adding the ifdef include statement', self.out_log, self.global_log)
154 f.write('\n')
155 f.write('; Include Position restraint file\n')
156 f.write('#ifdef '+str(posre_name)+'\n')
157 f.write('#include "'+str(Path(self.io_dict['out'].get("output_itp_path", "")).name)+'"\n')
158 f.write('#endif\n')
159 f.write('\n')
161 # For single-chain topologies
162 # Including new ITP in the TOP file
163 if not multi_chain:
165 # Read all lines of the top file
166 with open(top_file, 'r') as f:
167 lines = f.readlines()
169 main_chain = False
170 index = 0
172 # Find the index of the line where the custom position restraints are going to be included
173 for line in lines:
175 # Find the moleculetype directive of the main chain
176 if line.startswith('Protein_chain_'+chain):
177 main_chain = True
178 index = lines.index(line) + 3
180 # Find the end of the moleculetype directive of the main chain
181 if main_chain:
182 if line.startswith('[system]') or line.startswith('[molecules]') or line.startswith('#include ') or line.startswith('#ifdef POSRES'):
183 index = lines.index(line) - 1
184 break
186 if index == 0:
187 raise ValueError(f"Protein_chain_{chain} not found in the topology file")
189 # Include the custom position restraints in the top file
190 lines.insert(index, '\n')
191 lines.insert(index + 1, '; Include Position restraint file\n')
192 lines.insert(index + 2, '#ifdef '+str(posre_name)+'\n')
193 lines.insert(index + 3, '#include "'+str(Path(self.io_dict['out'].get("output_itp_path", "")).name)+'"\n')
194 lines.insert(index + 4, '#endif\n')
196 # Write the new top file
197 with open(top_file, 'w') as f:
198 f.writelines(lines)
200 # zip topology
201 fu.zip_top(zip_file=self.io_dict['out'].get("output_top_zip_path", ""), top_file=top_file, out_log=self.out_log, remove_original_files=self.remove_tmp)
203 # Remove temporal files
204 self.remove_tmp_files()
206 self.check_arguments(output_files_created=True, raise_exception=False)
207 return 0
210def ndx2resttop(input_ndx_path: str, input_top_zip_path: str, output_top_zip_path: str,
211 properties: Optional[dict] = None, **kwargs) -> int:
212 """Create :class:`Ndx2resttop <gromacs_extra.ndx2resttop.Ndx2resttop>` class and
213 execute the :meth:`launch() <gromacs_extra.ndx2resttop.Ndx2resttop.launch>` method."""
214 return Ndx2resttop(input_ndx_path=input_ndx_path,
215 input_top_zip_path=input_top_zip_path,
216 output_top_zip_path=output_top_zip_path,
217 properties=properties, **kwargs).launch()
220ndx2resttop.__doc__ = Ndx2resttop.__doc__
223def main():
224 """Command line execution of this building block. Please check the command line documentation."""
225 parser = argparse.ArgumentParser(description="Wrapper for the GROMACS extra ndx2resttop module.",
226 formatter_class=lambda prog: argparse.RawTextHelpFormatter(prog, width=99999))
227 parser.add_argument('-c', '--config', required=False, help="This file can be a YAML file, JSON file or JSON string")
229 # Specific args of each building block
230 required_args = parser.add_argument_group('required arguments')
231 required_args.add_argument('--input_ndx_path', required=True)
232 required_args.add_argument('--input_top_zip_path', required=True)
233 required_args.add_argument('--output_top_zip_path', required=True)
235 args = parser.parse_args()
236 config = args.config if args.config else None
237 properties = settings.ConfReader(config=config).get_prop_dic()
239 # Specific call of each building block
240 ndx2resttop(input_ndx_path=args.input_ndx_path, input_top_zip_path=args.input_top_zip_path,
241 output_top_zip_path=args.output_top_zip_path, properties=properties)
244if __name__ == '__main__':
245 main()