Author: bugman Date: Fri Oct 6 08:55:01 2006 New Revision: 2601 URL: http://svn.gna.org/viewcvs/relax?rev=2601&view=rev Log: Splitting of the 'pdb()' user function into 'pdb.read()' and 'pdb.vectors()' (task #3838).
These changes complete task #3838 (https://gna.org/task/?3838). The details of the changes were discussed at https://mail.gna.org/public/relax-users/2006-10/msg00002.html (Message-id: <[EMAIL PROTECTED]>). The only difference between the post proposing the changes and the actual changes is that the proposed user function 'pdb.set_vectors()' is actually now called 'pdb.vectors()'. Modified: 1.3/generic_fns/pdb.py 1.3/prompt/interpreter.py 1.3/prompt/pdb.py Modified: 1.3/generic_fns/pdb.py URL: http://svn.gna.org/viewcvs/relax/1.3/generic_fns/pdb.py?rev=2601&r1=2600&r2=2601&view=diff ============================================================================== --- 1.3/generic_fns/pdb.py (original) +++ 1.3/generic_fns/pdb.py Fri Oct 6 08:55:01 2006 @@ -114,7 +114,7 @@ i = i + 1 - def load(self, run=None, file=None, dir=None, model=None, heteronuc=None, proton=None, load_seq=1, calc_vectors=1, fail=1, print_flag=1): + def read(self, run=None, file=None, dir=None, model=None, load_seq=1, fail=1, print_flag=1): """The pdb loading function.""" # Arguments. @@ -122,10 +122,7 @@ self.file = file self.dir = dir self.model = model - self.heteronuc = heteronuc - self.proton = proton self.load_seq = load_seq - self.calc_vectors = calc_vectors self.fail = fail self.print_flag = print_flag @@ -193,10 +190,6 @@ # Load into Molmol (if running). self.relax.generic.molmol.open_pdb(self.run) - # Calculate the unit XH vectors. - if calc_vectors: - self.vectors() - def set_vector(self, run=None, res=None, xh_vect=None): """Function for setting the XH unit vectors.""" @@ -205,8 +198,34 @@ self.relax.data.res[run][res].xh_vect = xh_vect - def vectors(self): - """Function for calculating the XH unit vector from the loaded structure.""" + def vectors(self, heteronuc=None, proton=None, res_num=None, res_name=None): + """Function for calculating/extracting the XH unit vector from the loaded structure.""" + + # Arguments. + self.heteronuc = heteronuc + self.proton = proton + + # Test if the PDB file has been loaded. + if not self.relax.data.pdb.has_key(run): + raise RelaxPdbError, run + + # Test if sequence data is loaded. + if not self.relax.data.res.has_key(run): + raise RelaxNoSequenceError, run + + # Test if the residue number is a valid regular expression. + if type(num) == str: + try: + compile(num) + except: + raise RelaxRegExpError, ('residue number', num) + + # Test if the residue name is a valid regular expression. + if name: + try: + compile(name) + except: + raise RelaxRegExpError, ('residue name', name) # Print out. if self.print_flag: @@ -235,24 +254,44 @@ # Loop over the sequence. for j in xrange(len(self.relax.data.res[self.run])): + # Remap the data structure 'self.relax.data.res[self.run][j]'. + data = self.relax.data.res[self.run][j] + + # Skip unselected residues. + if not data.select: + continue + + # Skip the residue if there is no match to 'num'. + if type(num) == int: + if not data.num == num: + continue + elif type(num) == str: + if not match(num, `data.num`): + continue + + # Skip the residue if there is no match to 'name'. + if name != None: + if not match(name, data.name): + continue + # Find the corresponding residue in the PDB. pdb_res = None for k in xrange(len(pdb_residues)): - if self.relax.data.res[self.run][j].num == pdb_residues[k].number: + if data.num == pdb_residues[k].number: pdb_res = pdb_residues[k] break if pdb_res == None: - raise RelaxNoResError, self.relax.data.res[self.run][j].num + raise RelaxNoResError, data.num # Test if the proton atom exists for residue i. if not pdb_res.atoms.has_key(self.proton): - warn(RelaxNoAtomWarning(self.proton, self.relax.data.res[self.run][j].num)) - self.relax.data.res[self.run][j].xh_vect.append(None) + warn(RelaxNoAtomWarning(self.proton, data.num)) + data.xh_vect.append(None) # Test if the heteronucleus atom exists for residue i. elif not pdb_res.atoms.has_key(self.heteronuc): - warn(RelaxNoAtomWarning(self.heteronuc, self.relax.data.res[self.run][j].num)) - self.relax.data.res[self.run][j].xh_vect.append(None) + warn(RelaxNoAtomWarning(self.heteronuc, data.num)) + data.xh_vect.append(None) # Calculate the vector. else: @@ -271,12 +310,12 @@ # Test for zero length. if norm_factor == 0.0: if self.print_flag: - print "The XH bond vector for residue " + `self.relax.data.res[self.run][j].num` + " is of zero length." - self.relax.data.res[self.run][j].xh_vect.append(None) + print "The XH bond vector for residue " + `data.num` + " is of zero length." + data.xh_vect.append(None) # Calculate the normalised vector. else: - self.relax.data.res[self.run][j].xh_vect.append(vector / norm_factor) + data.xh_vect.append(vector / norm_factor) # Print out. if self.print_flag: @@ -287,9 +326,29 @@ # Average the vectors and convert xh_vect from an array of vectors to a vector. for i in xrange(len(self.relax.data.res[self.run])): + # Remap the data structure 'self.relax.data.res[self.run][j]'. + data = self.relax.data.res[self.run][j] + + # Skip unselected residues. + if not data.select: + continue + + # Skip the residue if there is no match to 'num'. + if type(num) == int: + if not data.num == num: + continue + elif type(num) == str: + if not match(num, `data.num`): + continue + + # Skip the residue if there is no match to 'name'. + if name != None: + if not match(name, data.name): + continue + # No vectors. - if self.relax.data.res[self.run][i].xh_vect[0] == None: - del self.relax.data.res[self.run][i].xh_vect + if data.xh_vect[0] == None: + del data.xh_vect continue # Average vectors. @@ -298,10 +357,10 @@ # Sum the vectors. for j in xrange(num_str): # Sum. - ave_vector = ave_vector + self.relax.data.res[self.run][i].xh_vect[j] + ave_vector = ave_vector + data.xh_vect[j] # Average the vector. ave_vector = ave_vector / num_str # Replace the temporary vector list with the normalised average vector. - self.relax.data.res[self.run][i].xh_vect = ave_vector / sqrt(dot(ave_vector, ave_vector)) + data.xh_vect = ave_vector / sqrt(dot(ave_vector, ave_vector)) Modified: 1.3/prompt/interpreter.py URL: http://svn.gna.org/viewcvs/relax/1.3/prompt/interpreter.py?rev=2601&r1=2600&r2=2601&view=diff ============================================================================== --- 1.3/prompt/interpreter.py (original) +++ 1.3/prompt/interpreter.py Fri Oct 6 08:55:01 2006 @@ -47,7 +47,6 @@ from minimisation import Minimisation from model_selection import Modsel from nuclei import Nuclei -from pdb import PDB # User classes. from dasha import Dasha @@ -59,6 +58,7 @@ from monte_carlo import Monte_carlo from noe import Noe from palmer import Palmer +from pdb import PDB from relax_data import Relax_data from relax_fit import Relax_fit from results import Results @@ -102,7 +102,6 @@ self._Modsel = Modsel(relax) self._Nuclei = Nuclei(relax) self._OpenDX = OpenDX(relax) - self._PDB = PDB(relax) self._system = system # Place the user classes into the interpreter class namespace. @@ -116,6 +115,7 @@ self._Monte_carlo = Monte_carlo(relax) self._Noe = Noe(relax) self._Palmer = Palmer(relax) + self._PDB = PDB(relax) self._Relax_data = Relax_data(relax) self._Relax_fit = Relax_fit(relax) self._Results = Results(relax) @@ -160,7 +160,6 @@ minimise = self._Minimisation.minimise model_selection = self._Modsel.model_selection nuclei = self._Nuclei.nuclei - pdb = self._PDB.pdb # Place the user classes in the local namespace. dasha = self._Dasha @@ -173,6 +172,7 @@ monte_carlo = self._Monte_carlo noe = self._Noe palmer = self._Palmer + pdb = self._PDB relax_data = self._Relax_data relax_fit = self._Relax_fit results = self._Results Modified: 1.3/prompt/pdb.py URL: http://svn.gna.org/viewcvs/relax/1.3/prompt/pdb.py?rev=2601&r1=2600&r2=2601&view=diff ============================================================================== --- 1.3/prompt/pdb.py (original) +++ 1.3/prompt/pdb.py Fri Oct 6 08:55:01 2006 @@ -1,6 +1,6 @@ ############################################################################### # # -# Copyright (C) 2003, 2004 Edward d'Auvergne # +# Copyright (C) 2003, 2004, 2006 Edward d'Auvergne # # # # This file is part of the program relax. # # # @@ -22,15 +22,23 @@ import sys +import help + class PDB: def __init__(self, relax): - """Class containing the function for loading a pdb file.""" - - self.relax = relax - - - def pdb(self, run=None, file=None, dir=None, model=None, heteronuc='N', proton='H', load_seq=1): + # Help. + self.__relax_help__ = \ + """Class containing the PDB related functions.""" + + # Add the generic help string. + self.__relax_help__ = self.__relax_help__ + "\n" + help.relax_class_help + + # Place relax in the class namespace. + self.__relax__ = relax + + + def read(self, run=None, file=None, dir=None, model=None, heteronuc='N', proton='H', load_seq=1): """The pdb loading function. Keyword Arguments @@ -43,10 +51,6 @@ dir: The directory where the file is located. model: The PDB model number. - - heteronuc: The name of the heteronucleus as specified in the PDB file. - - proton: The name of the proton as specified in the PDB file. load_seq: A flag specifying whether the sequence should be loaded from the PDB file. @@ -61,39 +65,31 @@ To load the sequence from the PDB file, set the 'load_seq' flag to 1. If the sequence has previously been loaded, then this flag will be ignored. - Once the PDB structures are loaded, unit XH bond vectors will be calculated. The vectors - are calculated using the atomic coordinates of the atoms specified by the arguments - heteronuc and proton. If more than one model structure is loaded, the unit XH vectors for - each model will be calculated and the final unit XH vector will be taken as the average. - Example ~~~~~~~ To load all structures from the PDB file 'test.pdb' in the directory '~/pdb' for use in the - model-free analysis run 'm8' where the heteronucleus in the PDB file is 'N' and the proton - is 'H', type: - - relax> pdb('m8', 'test.pdb', '~/pdb', 1, 'N', 'H') - relax> pdb(run='m8', file='test.pdb', dir='pdb', model=1, heteronuc='N', proton='H') + model-free analysis run 'm8', type: + + relax> pdb.read('m8', 'test.pdb', '~/pdb', 1) + relax> pdb.read(run='m8', file='test.pdb', dir='pdb', model=1) To load the 10th model from the file 'test.pdb', use: - relax> pdb('m1', 'test.pdb', model=10) - relax> pdb(run='m1', file='test.pdb', model=10) + relax> pdb.read('m1', 'test.pdb', model=10) + relax> pdb.read(run='m1', file='test.pdb', model=10) """ # Function intro text. if self.relax.interpreter.intro: - text = sys.ps3 + "pdb(" + text = sys.ps3 + "pdb.read(" text = text + "run=" + `run` text = text + ", file=" + `file` text = text + ", dir=" + `dir` text = text + ", model=" + `model` - text = text + ", heteronuc=" + `heteronuc` - text = text + ", proton=" + `proton` text = text + ", load_seq=" + `load_seq` + ")" print text @@ -113,6 +109,79 @@ if model != None and type(model) != int: raise RelaxIntError, ('model', model) + # The load sequence argument. + if type(load_seq) != int or (load_seq != 0 and load_seq != 1): + raise RelaxBinError, ('load sequence flag', load_seq) + + # Execute the functional code. + self.__relax__.generic.pdb.read(run=run, file=file, dir=dir, model=model, load_seq=load_seq) + + + def vectors(self, run=None, heteronuc='N', proton='H', res_num=None, res_name=None): + """Function for calculating/extracting XH vectors from the structure. + + Keyword arguments + ~~~~~~~~~~~~~~~~~ + + run: The run to assign the vectors to. + + heteronuc: The heteronucleus name as specified in the PDB file. + + proton: The name of the proton as specified in the PDB file. + + res_num: The residue number. + + res_name: The name of the residue. + + + Description + ~~~~~~~~~~~ + + Once the PDB structures have been loaded, the unit XH bond vectors must be calculated for + the non-spherical diffusion models. The vectors are calculated using the atomic coordinates + of the atoms specified by the arguments heteronuc and proton. If more than one model + structure is loaded, the unit XH vectors for each model will be calculated and the final + unit XH vector will be taken as the average. + + + Example + ~~~~~~~ + + To calculate the XH vectors of the backbone amide nitrogens where in the PDB file the + backbone nitrogen is called 'N' and the attached proton is called 'H', assuming the run + 'test', type: + + relax> pdb.vectors('test') + relax> pdb.vectors('test', 'N') + relax> pdb.vectors('test', 'N', 'H') + relax> pdb.vectors('test', heteronuc='N', proton='H') + + If the attached proton is called 'HN', type: + + relax> pdb.vectors('test', proton='HN') + + If you are working with RNA, you can use the residue name identifier to calculate the + vectors for each residue separately. For example: + + relax> pdb.vectors('m1', 'N1', 'H1', res_name='G') + relax> pdb.vectors('m1', 'N3', 'H3', res_name='U') + + """ + + # Function intro text. + if self.relax.interpreter.intro: + text = sys.ps3 + "pdb.vectors(" + text = text + "run=" + `run` + text = text + ", heteronuc=" + `heteronuc` + text = text + ", proton=" + `proton` + text = text + ", res_num=" + `res_num` + text = text + ", res_name=" + `res_name` + ")" + print text + + # The run argument. + if type(run) != str: + raise RelaxStrError, ('run', run) + # The heteronucleus argument. if type(heteronuc) != str: raise RelaxStrError, ('heteronucleus', heteronuc) @@ -121,9 +190,13 @@ if type(proton) != str: raise RelaxStrError, ('proton', proton) - # The load sequence argument. - if type(load_seq) != int or (load_seq != 0 and load_seq != 1): - raise RelaxBinError, ('load sequence flag', load_seq) + # Residue number. + if type(res_num) != int: + raise RelaxIntError, ('residue number', res_num) + + # Residue name. + if type(res_name) != str: + raise RelaxStrError, ('residue name', res_name) # Execute the functional code. - self.relax.generic.pdb.load(run=run, file=file, dir=dir, model=model, heteronuc=heteronuc, proton=proton, load_seq=load_seq) + self.__relax__.generic.pdb.vectors(run=run, heteronuc=heteronuc, proton=proton, res_num=res_num, res_name=res_name) _______________________________________________ relax (http://nmr-relax.com) This is the relax-commits mailing list relax-commits@gna.org To unsubscribe from this list, get a password reminder, or change your subscription options, visit the list information page at https://mail.gna.org/listinfo/relax-commits