Package generic_fns :: Module dasha
[hide private]
[frames] | no frames]

Source Code for Module generic_fns.dasha

  1  ############################################################################### 
  2  #                                                                             # 
  3  # Copyright (C) 2005 Edward d'Auvergne                                        # 
  4  #                                                                             # 
  5  # This file is part of the program relax.                                     # 
  6  #                                                                             # 
  7  # relax is free software; you can redistribute it and/or modify               # 
  8  # it under the terms of the GNU General Public License as published by        # 
  9  # the Free Software Foundation; either version 2 of the License, or           # 
 10  # (at your option) any later version.                                         # 
 11  #                                                                             # 
 12  # relax is distributed in the hope that it will be useful,                    # 
 13  # but WITHOUT ANY WARRANTY; without even the implied warranty of              # 
 14  # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the               # 
 15  # GNU General Public License for more details.                                # 
 16  #                                                                             # 
 17  # You should have received a copy of the GNU General Public License           # 
 18  # along with relax; if not, write to the Free Software                        # 
 19  # Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA   # 
 20  #                                                                             # 
 21  ############################################################################### 
 22   
 23  from math import pi 
 24  from os import F_OK, access, chdir, getcwd, system 
 25  from re import match, search 
 26  from string import split 
 27  import sys 
 28   
 29   
30 -class Dasha:
31 - def __init__(self, relax):
32 """Class used to create and process input and output for the program Modelfree 4.""" 33 34 self.relax = relax
35 36
37 - def create(self, run=None, algor='LM', dir=None, force=0):
38 """Function for creating the Dasha script file 'dir/dasha_script'.""" 39 40 # Arguments. 41 self.run = run 42 self.algor = algor 43 self.dir = dir 44 self.force = force 45 46 # Test if the run exists. 47 if not self.run in self.relax.data.run_names: 48 raise RelaxNoRunError, self.run 49 50 # Test if sequence data is loaded. 51 if not self.relax.data.res.has_key(self.run): 52 raise RelaxNoSequenceError, self.run 53 54 # Determine the parameter set. 55 self.param_set = self.relax.specific.model_free.determine_param_set_type(self.run) 56 57 # Test if diffusion tensor data for the run exists. 58 if self.param_set != 'local_tm' and not self.relax.data.diff.has_key(self.run): 59 raise RelaxNoTensorError, self.run 60 61 # Test if the PDB file has been loaded (for the spheroid and ellipsoid). 62 if self.param_set != 'local_tm' and self.relax.data.diff[self.run].type != 'sphere' and not self.relax.data.pdb.has_key(self.run): 63 raise RelaxNoPdbError, self.run 64 65 # Test if the nucleus type has been set. 66 if not hasattr(self.relax.data, 'gx'): 67 raise RelaxNucleusError 68 69 # Test the optimisation algorithm. 70 if algor not in ['LM', 'NR']: 71 raise RelaxError, "The Dasha optimisation algorithm " + `algor` + " is unknown, it should either be 'LM' or 'NR'." 72 73 # Directory creation. 74 if self.dir == None: 75 self.dir = self.run 76 self.relax.IO.mkdir(self.dir, print_flag=0) 77 78 # Number of field strengths and values. 79 self.num_frq = 0 80 self.frq = [] 81 for i in xrange(len(self.relax.data.res[self.run])): 82 if hasattr(self.relax.data.res[self.run][i], 'num_frq'): 83 if self.relax.data.res[self.run][i].num_frq > self.num_frq: 84 # Number of field strengths. 85 self.num_frq = self.relax.data.res[self.run][i].num_frq 86 87 # Field strength values. 88 for frq in self.relax.data.res[self.run][i].frq: 89 if frq not in self.frq: 90 self.frq.append(frq) 91 92 # Calculate the angle alpha of the XH vector in the spheroid diffusion frame. 93 if self.relax.data.diff[self.run].type == 'spheroid': 94 self.relax.generic.angles.spheroid_frame(self.run) 95 96 # Calculate the angles theta and phi of the XH vector in the ellipsoid diffusion frame. 97 elif self.relax.data.diff[self.run].type == 'ellipsoid': 98 self.relax.generic.angles.ellipsoid_frame(self.run) 99 100 # The 'dasha_script' file. 101 script = self.relax.IO.open_write_file(file_name='dasha_script', dir=self.dir, force=self.force) 102 self.create_script(script) 103 script.close()
104 105
106 - def create_script(self, file):
107 """Create the Dasha script file.""" 108 109 # Delete all data. 110 file.write('# Delete all data.\n') 111 file.write('del 1 10000\n') 112 113 # Nucleus type. 114 file.write('\n# Nucleus type.\n') 115 nucleus = self.relax.generic.nuclei.find_nucleus() 116 if nucleus == 'N': 117 nucleus = 'N15' 118 elif nucleus == 'C': 119 nucleus = 'C13' 120 else: 121 raise RelaxError, 'Cannot handle the nucleus type ' + `nucleus` + ' within Dasha.' 122 file.write('set nucl ' + nucleus + '\n') 123 124 # Number of frequencies. 125 file.write('\n# Number of frequencies.\n') 126 file.write('set n_freq ' + `self.relax.data.num_frq[self.run]` + '\n') 127 128 # Frequency values. 129 file.write('\n# Frequency values.\n') 130 for i in xrange(self.relax.data.num_frq[self.run]): 131 file.write('set H1_freq ' + `self.relax.data.frq[self.run][i] / 1e6` + ' ' + `i+1` + '\n') 132 133 # Set the diffusion tensor. 134 file.write('\n# Set the diffusion tensor.\n') 135 if self.param_set != 'local_tm': 136 # Sphere. 137 if self.relax.data.diff[self.run].type == 'sphere': 138 file.write('set tr ' + `self.relax.data.diff[self.run].tm / 1e-9` + '\n') 139 140 # Spheroid. 141 elif self.relax.data.diff[self.run].type == 'spheroid': 142 file.write('set tr ' + `self.relax.data.diff[self.run].tm / 1e-9` + '\n') 143 144 # Ellipsoid. 145 elif self.relax.data.diff[self.run].type == 'ellipsoid': 146 # Get the eigenvales. 147 Dx, Dy, Dz = self.relax.generic.diffusion_tensor.return_eigenvalues(self.run) 148 149 # Geometric parameters. 150 file.write('set tr ' + `self.relax.data.diff[self.run].tm / 1e-9` + '\n') 151 file.write('set D1/D3 ' + `Dx / Dz` + '\n') 152 file.write('set D2/D3 ' + `Dy / Dz` + '\n') 153 154 # Orientational parameters. 155 file.write('set alfa ' + `self.relax.data.diff[self.run].alpha / (2.0 * pi) * 360.0` + '\n') 156 file.write('set betta ' + `self.relax.data.diff[self.run].beta / (2.0 * pi) * 360.0` + '\n') 157 file.write('set gamma ' + `self.relax.data.diff[self.run].gamma / (2.0 * pi) * 360.0` + '\n') 158 159 # Reading the relaxation data. 160 file.write('\n# Reading the relaxation data.\n') 161 file.write('echo Reading the relaxation data.\n') 162 noe_index = 1 163 r1_index = 1 164 r2_index = 1 165 for i in xrange(self.relax.data.num_ri[self.run]): 166 # NOE. 167 if self.relax.data.ri_labels[self.run][i] == 'NOE': 168 # Data set number. 169 number = noe_index 170 171 # Data type. 172 data_type = 'noe' 173 174 # Increment the data set index. 175 noe_index = noe_index + 1 176 177 # R1. 178 elif self.relax.data.ri_labels[self.run][i] == 'R1': 179 # Data set number. 180 number = r1_index 181 182 # Data type. 183 data_type = '1/T1' 184 185 # Increment the data set index. 186 r1_index = r1_index + 1 187 188 # R2. 189 elif self.relax.data.ri_labels[self.run][i] == 'R2': 190 # Data set number. 191 number = r2_index 192 193 # Data type. 194 data_type = '1/T2' 195 196 # Increment the data set index. 197 r2_index = r2_index + 1 198 199 # Set the data type. 200 if number == 1: 201 file.write('\nread < ' + data_type + '\n') 202 else: 203 file.write('\nread < ' + data_type + ' ' + `number` + '\n') 204 205 # The relaxation data. 206 for j in xrange(len(self.relax.data.res[self.run])): 207 # Reassign the data. 208 data = self.relax.data.res[self.run][j] 209 210 # Data and errors. 211 file.write(`data.num` + ' ' + `data.relax_data[i]` + ' ' + `data.relax_error[i]` + '\n') 212 213 # Terminate the reading. 214 file.write('exit\n') 215 216 # Individual residue optimisation. 217 if self.param_set == 'mf': 218 # Loop over the residues. 219 for i in xrange(len(self.relax.data.res[self.run])): 220 # Reassign the data. 221 data = self.relax.data.res[self.run][i] 222 223 # Comment. 224 file.write('\n\n\n# Residue ' + `data.num` + '\n\n') 225 226 # Echo. 227 file.write('echo Optimisation of residue ' + `data.num` + '\n') 228 229 # Select the residue. 230 file.write('\n# Select the residue.\n') 231 file.write('set cres ' + `data.num` + '\n') 232 233 # The angle alpha of the XH vector in the spheroid diffusion frame. 234 if self.relax.data.diff[self.run].type == 'spheroid': 235 file.write('set teta ' + `data.alpha` + '\n') 236 237 # The angles theta and phi of the XH vector in the ellipsoid diffusion frame. 238 elif self.relax.data.diff[self.run].type == 'ellipsoid': 239 file.write('\n# Setting the spherical angles of the XH vector in the ellipsoid diffusion frame.\n') 240 file.write('set teta ' + `data.theta` + '\n') 241 file.write('set fi ' + `data.phi` + '\n') 242 243 # The 'jmode'. 244 if 'ts' in data.params: 245 jmode = 3 246 elif 'te' in data.params: 247 jmode = 2 248 elif 'S2' in data.params: 249 jmode = 1 250 251 # Chemical exchange. 252 if 'Rex' in data.params: 253 exch = 1 254 else: 255 exch = 0 256 257 # Anisotropic diffusion. 258 if self.relax.data.diff[self.run].type == 'sphere': 259 anis = 0 260 else: 261 anis = 1 262 263 # Axial symmetry. 264 if self.relax.data.diff[self.run].type == 'spheroid': 265 sym = 1 266 else: 267 sym = 0 268 269 # Set the jmode. 270 file.write('\n# Set the jmode.\n') 271 file.write('set def jmode ' + `jmode`) 272 if exch: 273 file.write(' exch') 274 if anis: 275 file.write(' anis') 276 if sym: 277 file.write(' sym') 278 file.write('\n') 279 280 # Fix the tf parameter if it isn't in the model. 281 if not 'tf' in data.params and jmode == 3: 282 file.write('\n# Fix the tf parameter.\n') 283 file.write('fix tf 0\n') 284 285 # Bond length. 286 file.write('\n# Bond length.\n') 287 file.write('set r_hx ' + `data.r / 1e-10` + '\n') 288 289 # CSA value. 290 file.write('\n# CSA value.\n') 291 file.write('set csa ' + `data.csa / 1e-6` + '\n') 292 293 # Parameter default values. 294 file.write('\n# Parameter default values.\n') 295 file.write('reset jmode ' + `data.num` + '\n') 296 297 # Optimisation of all residues. 298 file.write('\n\n\n# Optimisation of all residues.\n') 299 if self.algor == 'LM': 300 file.write('lmin all') 301 elif self.algor == 'NR': 302 file.write('min all') 303 304 # Show the results. 305 file.write('\n# Show the results.\n') 306 file.write('echo\n') 307 file.write('show all\n') 308 309 # Write the results. 310 file.write('\n# Write the results.\n') 311 file.write('write S2.out S\n') 312 file.write('write S2f.out Sf\n') 313 file.write('write S2s.out Ss\n') 314 file.write('write te.out te\n') 315 file.write('write tf.out tf\n') 316 file.write('write ts.out ts\n') 317 file.write('write Rex.out rex\n') 318 file.write('write chi2.out F\n') 319 320 else: 321 raise RelaxError, 'Optimisation of the parameter set ' + `self.param_set` + ' currently not supported.'
322 323
324 - def execute(self, run, dir, force):
325 """Function for executing Dasha.""" 326 327 # The current directory. 328 orig_dir = getcwd() 329 330 # The directory. 331 if dir == None: 332 dir = run 333 if not access(dir, F_OK): 334 raise RelaxDirError, ('Dasha', dir) 335 336 # Change to this directory. 337 chdir(dir) 338 339 # Catch failures and return to the correct directory. 340 try: 341 # Test if the 'dasha_script' script file exists. 342 if not access('dasha_script', F_OK): 343 raise RelaxFileError, ('dasha script', 'dasha_script') 344 345 # Execute Dasha. 346 system('dasha < dasha_script | tee dasha_results') 347 348 # Failure. 349 except: 350 # Change back to the original directory. 351 chdir(orig_dir) 352 353 # Reraise the error. 354 raise 355 356 # Change back to the original directory. 357 chdir(orig_dir) 358 359 # Print some blank lines (aesthetics) 360 sys.stdout.write('\n\n')
361 362
363 - def extract(self, run, dir):
364 """Function for extracting the Dasha results out of the 'dasha_results' file.""" 365 366 # Arguments. 367 self.run = run 368 369 # Test if sequence data is loaded. 370 if not self.relax.data.res.has_key(self.run): 371 raise RelaxNoSequenceError, self.run 372 373 # The directory. 374 if dir == None: 375 dir = run 376 if not access(dir, F_OK): 377 raise RelaxDirError, ('Dasha', dir) 378 379 # Loop over the parameters. 380 for param in ['S2', 'S2f', 'S2s', 'te', 'tf', 'ts', 'Rex']: 381 # The file name. 382 file_name = dir + '/' + param + '.out' 383 384 # Test if the file exists. 385 if not access(file_name, F_OK): 386 raise RelaxFileError, ('Dasha', file_name) 387 388 # Scaling. 389 if param in ['te', 'tf', 'ts']: 390 scaling = 1e-9 391 elif param == 'Rex': 392 scaling = 1.0 / (2.0 * pi * self.relax.data.frq[self.run][0]) ** 2 393 else: 394 scaling = 1.0 395 396 # Set the values. 397 self.relax.generic.value.read(self.run, param=param, scaling=scaling, file=file_name, num_col=0, name_col=None, data_col=1, error_col=2) 398 399 # Extract the chi-squared values. 400 file_name = dir + '/chi2.out' 401 402 # Test if the file exists. 403 if not access(file_name, F_OK): 404 raise RelaxFileError, ('Dasha', file_name) 405 406 # Set the values. 407 self.relax.generic.value.read(self.run, param='chi2', file=file_name, num_col=0, name_col=None, data_col=1, error_col=2)
408