Source Code for Module Biskit.AmberRstParser

  1  ## Automatically adapted for numpy.oldnumeric Mar 26, 2007 by alter_code1.py 
  2   
  3  ## 
  4  ## Biskit, a toolkit for the manipulation of macromolecular structures 
  5  ## Copyright (C) 2004-2006 Raik Gruenberg & Johan Leckner 
  6  ## 
  7  ## This program is free software; you can redistribute it and/or 
  8  ## modify it under the terms of the GNU General Public License as 
  9  ## published by the Free Software Foundation; either version 2 of the 
 10  ## License, or any later version. 
 11  ## 
 12  ## This program 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 GNU 
 15  ## General Public License for more details. 
 16  ## 
 17  ## You find a copy of the GNU General Public License in the file 
 18  ## license.txt along with this program; if not, write to the Free 
 19  ## Software Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. 
 20  ## 
 21  ## 
 22  ## $Revision: 2.7 $ 
 23  ## last $Date: 2007/03/26 18:40:40 $ 
 24  ## last $Author: graik $ 
 25  """ 
 26  Parse Amber restart files. 
 27  """ 
 28       
 29  import re 
 30  import numpy.oldnumeric as N 
 31  import os.path 
 32   
 33  from AmberCrdParser import ParseError 
 34  from PDBModel import PDBModel 
 35  import tools as T 
 36   
 37 -class AmberRstParser: 
 38      """Convert an Amber restart file to array, PDBModel or a Amber crd file. 
 39   
 40      Note: AmberRstParser is currently ignoring both the velocity and 
 41      boxinfo record (although this could be easily changed). 
 42      """ 
 43   
 44 -    def __init__( self, frst ): 
 45          """ 
 46          @param frst: input restart file 
 47          @type  frst: str 
 48          """ 
 49          self.frst = T.absfile( frst ) 
 50          self.crd  = open( self.frst ) 
 51   
 52          self.n = 0       #: number of atoms 
 53          self.lines_per_frame = 0 
 54          self.xyz = None  #: will hold coordinate array 
 55          self.box = None  #: will hold box array if any 
 56   
 57          ## pre-compile pattern for line2numbers 
 58          xnumber = "-*\d+\.\d+"              # optionally negtive number 
 59          xspace  = ' *'                      # one or more space char 
 60          self.xnumbers = re.compile('('+xspace+xnumber+')') 
 61   
 62   
 63 -    def __del__(self): 
 64          try: 
 65              self.crd.close() 
 66          except: 
 67              pass 
 68   
 69   
 70 -    def __nextLine( self ): 
 71          """Extract next line of coordinates from crd file 
 72   
 73          @return: coordinates 
 74          @rtype: [float]     
 75          """ 
 76          l = self.crd.readline() 
 77          if l == '': 
 78              raise EOFError('EOF') 
 79   
 80          match = self.xnumbers.findall( l ) 
 81          return [ round( float(strCrd),7) for strCrd in match ]  
 82   
 83   
 84 -    def __frame( self ): 
 85          """Collect next complete coordinate frame 
 86   
 87          @return: coordinate frame 
 88          @rtype: array 
 89          """ 
 90          self.xyz = [ self.__nextLine() for i in range(self.lines_per_frame) ] 
 91   
 92          return N.reshape(self.xyz, ( self.n, 3 ) ).astype(N.Float32) 
 93   
 94   
 95 -    def getXyz( self ): 
 96          """Get coordinate array. 
 97           
 98          @return: coordinates, N.array( N x 3, 'f') 
 99          @rtype: array 
100           
101          @raise ParseError: if can't interprete second line 
102          """ 
103          if not self.xyz: 
104   
105              ## skip first empty line 
106              self.crd.readline() 
107   
108              try: 
109                  self.n, self.time = self.crd.readline().split() 
110                  self.n = int( self.n ) 
111                  self.time = float( self.time ) 
112              except: 
113                  raise ParseError("Can't interprete second line of "+self.frst) 
114   
115              ## pre-compute lines expected per frame 
116              self.lines_per_frame = self.n / 2 
117              if self.n % 2 != 0: 
118                  self.lines_per_frame += 1 
119   
120              self.xyz = self.__frame() 
121   
122          return self.xyz 
123   
124   
125 -    def getModel( self, ref, rnAmber=0 ): 
126          """ 
127          Get model. 
128           
129          @param ref: reference with same number and order of atoms 
130          @type  ref: PDBModel 
131          @param rnAmber: rename Amber to standard residues (HIE, HID, HIP, CYX) 
132          @type  rnAmber: 1|0 
133           
134          @return: PDBModel 
135          @rtype: PDBModel 
136          """ 
137          if not self.xyz: 
138              self.getXyz() 
139   
140          result = ref.clone() 
141          result.setXyz( self.xyz ) 
142          if rnAmber: 
143              result.renameAmberRes() 
144   
145          return result 
146   
147   
148 -    def getFirstCrdLine( self ): 
149          """ 
150          Return the first line of Amber crd. 
151           
152          @return: first line of Amber crd formatted coordinate block 
153          @rtype: str 
154          """ 
155          if not self.xyz: 
156              self.getXyz() 
157   
158          result = "" 
159          for x in N.ravel( self.xyz )[:10]: 
160              result += "%8.3f" % x 
161   
162          return result + "\n" 
163   
164   
165 -    def writeCrd( self, fcrd, append=1, lastAtom=None ): 
166          """ 
167          Write/Append Amber-formatted block of coordinates to a file. 
168          If a file handle is given, the file will not be closed. 
169           
170          @param fcrd: file to write to 
171          @type  fcrd: str or file object 
172          @param append: append to existing file (default: 1) 
173          @type  append: 0|1 
174          @param lastAtom: skip all atoms beyond this one (default: None) 
175          @type  lastAtom: int 
176          """ 
177          if not self.xyz: 
178              self.getXyz() 
179   
180          if type( fcrd ) == file: 
181              ## take file handle 
182              f = fcrd 
183          else: 
184              ## create new file handle 
185              mode = 'w' 
186              if append: 
187                  mode = 'a' 
188              f = open( T.absfile( fcrd ), mode ) 
189               
190              newf = (mode=='w' or not os.path.exists( T.absfile(fcrd) )) 
191              if newf: 
192                  f.write("\n") 
193   
194          i = 0 
195          for x in N.ravel( self.xyz ): 
196              i = i + 1 
197   
198              f.write( "%8.3f" % x ) 
199   
200              if (i % 10) == 0: 
201                  f.write("\n") 
202   
203              if lastAtom and i / 3.0 == lastAtom: 
204                  break 
205   
206          if ((i % 10) != 0): 
207              f.write("\n") 
208   
209          if type( fcrd ) != file: 
210              ## don't close file that was already given 
211              f.close() 
212   
213  ###################### 
214  ### Module testing ### 
215  import Biskit.test as BT 
216   
217 -class Test(BT.BiskitTest): 
218      """Test AmberRstParser""" 
219   
220 -    def prepare(self): 
221          self.f    = T.testRoot()+'/amber/sim.rst' 
222          self.fref = T.testRoot()+'/amber/1HPT_0.pdb' 
223           
224          self.p = AmberRstParser( self.f ) 
225   
226 -    def test_getXyz( self ): 
227          """AmberRstParser.getXyz test""" 
228          self.xyz = self.p.getXyz() 
229          self.assertEqual( N.shape(self.xyz), (11200,3) ) 
230   
231 -class TestLong(Test): 
232      """long AmberRstParser test""" 
233      TAGS = [BT.LONG] 
234   
235 -    def test_getModel(self): 
236          """AmberRstParser.getModel test""" 
237          self.ref = PDBModel( self.fref ) 
238          self.model = self.p.getModel( self.ref ) 
239          self.assertEqual( len(self.model), 11200 ) 
240   
241   
242  if __name__ == '__main__': 
243   
244      ## run Test and push self.* fields into global namespace 
245      BT.localTest( ) 
246