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

1#!/usr/bin/env python3 

2 

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 

13 

14 

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. 

20 

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. 

29 

30 Examples: 

31 This is a use example of how to use the building block from Python:: 

32 

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) 

39 

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 """ 

48 

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 {} 

52 

53 # Call parent class constructor 

54 super().__init__(properties) 

55 self.locals_var_dict = locals().copy() 

56 

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 } 

62 

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') 

67 

68 # Check the properties 

69 self.check_properties(properties) 

70 self.check_arguments() 

71 

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 

78 

79 top_file = fu.unzip_top(zip_file=self.io_dict['in'].get("input_top_zip_path", ""), out_log=self.out_log) 

80 

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): 

85 

86 # New group 

87 if line.startswith('['): 

88 index_dic[line] = [index, 0] 

89 

90 # Update current group 

91 label = line 

92 

93 # Close previous group 

94 if index > 0: 

95 index_dic[label] = [index_dic[label][0], index] 

96 

97 # Last group of the file 

98 if index == len(lines)-1: 

99 index_dic[label] = [index_dic[label][0], index] 

100 

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 

105 

106 fu.log('Index_dic: '+str(index_dic), self.out_log, self.global_log) 

107 

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) 

110 

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) 

117 

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") 

121 

122 for triplet, posre_name in zip(self.ref_rest_chain_triplet_list, self.posres_names): 

123 

124 reference_group, restrain_group, chain = triplet 

125 

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') 

130 

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') 

144 

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') 

160 

161 # For single-chain topologies 

162 # Including new ITP in the TOP file 

163 if not multi_chain: 

164 

165 # Read all lines of the top file 

166 with open(top_file, 'r') as f: 

167 lines = f.readlines() 

168 

169 main_chain = False 

170 index = 0 

171 

172 # Find the index of the line where the custom position restraints are going to be included 

173 for line in lines: 

174 

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 

179 

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 

185 

186 if index == 0: 

187 raise ValueError(f"Protein_chain_{chain} not found in the topology file") 

188 

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') 

195 

196 # Write the new top file 

197 with open(top_file, 'w') as f: 

198 f.writelines(lines) 

199 

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) 

202 

203 # Remove temporal files 

204 self.remove_tmp_files() 

205 

206 self.check_arguments(output_files_created=True, raise_exception=False) 

207 return 0 

208 

209 

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() 

218 

219 

220ndx2resttop.__doc__ = Ndx2resttop.__doc__ 

221 

222 

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") 

228 

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) 

234 

235 args = parser.parse_args() 

236 config = args.config if args.config else None 

237 properties = settings.ConfReader(config=config).get_prop_dic() 

238 

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) 

242 

243 

244if __name__ == '__main__': 

245 main()