source: trunk/imports/G2pwd_fxye.py @ 1713

Last change on this file since 1713 was 1713, checked in by vondreele, 7 years ago

fix problems with imports of old TOF data; bank nos. now are as given in parm/data file rather than starting with 1...

  • Property svn:eol-style set to native
  • Property svn:keywords set to Date Author Revision URL Id
File size: 16.2 KB
Line 
1# -*- coding: utf-8 -*-
2########### SVN repository information ###################
3# $Date: 2015-03-16 16:03:41 +0000 (Mon, 16 Mar 2015) $
4# $Author: vondreele $
5# $Revision: 1713 $
6# $URL: trunk/imports/G2pwd_fxye.py $
7# $Id: G2pwd_fxye.py 1713 2015-03-16 16:03:41Z vondreele $
8########### SVN repository information ###################
9'''
10*Module G2pwd_fxye: GSAS data files*
11------------------------------------
12Routine to read in powder data in a variety of formats
13that are defined for GSAS.
14
15'''
16import sys
17import os.path as ospath
18import numpy as np
19import GSASIIIO as G2IO
20import GSASIIpath
21GSASIIpath.SetVersionNumber("$Revision: 1713 $")
22
23class GSAS_ReaderClass(G2IO.ImportPowderData):
24    'Routines to import powder data from a GSAS files'
25    def __init__(self):
26        super(self.__class__,self).__init__( # fancy way to self-reference
27            extensionlist=('.fxye','.raw','.gsas','.gda','gsa.','.RAW','.GSAS','.GDA','.GSA'),
28            strictExtension=False,
29            formatName = 'GSAS powder data',
30            longFormatName = 'GSAS powder data files (.fxye, .raw, .gsas...)'
31            )
32        self.clockWd = None
33
34    # Validate the contents -- look for a bank line
35    def ContentsValidator(self, filepointer):
36        'Validate by checking to see if the file has BANK lines'
37        #print 'ContentsValidator: '+self.formatName
38        for i,line in enumerate(filepointer):
39            self.GSAS = True
40            if i==0: # first line is always a comment
41                continue
42            if i==1 and line[:4].lower() == 'inst':
43                # 2nd line is optional instrument parameter file
44                continue
45            if line[0] == '#': continue
46            if line[:4] == 'BANK':
47                return True
48            elif line[:7] == 'Monitor': continue
49            elif line [:8] == 'TIME_MAP':          #LANSCE TOF data
50                return True
51            else:
52                self.errors = 'Unexpected information in line: '+str(i+1)
53                self.errors += '  '+str(line)
54                return False
55        return False # no bank records
56
57    def Reader(self,filename,filepointer, ParentFrame=None, **kwarg):
58        '''Read a GSAS (old formats) file of type FXY, FXYE, ESD or STD types.
59        If multiple datasets are requested, use self.repeat and buffer caching.
60        '''
61        def GetFXYEdata(File,Pos,Bank):
62            File.seek(Pos)
63            x = []
64            y = []
65            w = []
66            S = File.readline()
67            while S and S[:4] != 'BANK' and S[0] != '#':
68                vals = S.split()
69                x.append(float(vals[0])/100.)               #CW: from centidegrees to degrees
70                f = float(vals[1])
71                s = float(vals[2])
72                if f <= 0.0 or s <= 0.0:
73                    y.append(0.0)
74                    w.append(0.0)
75                else:
76                    y.append(float(vals[1]))
77                    w.append(1.0/float(vals[2])**2)
78                S = File.readline()
79            N = len(x)
80            return [np.array(x),np.array(y),np.array(w),np.zeros(N),np.zeros(N),np.zeros(N)]   
81           
82        def GetFXYdata(File,Pos,Bank):
83            File.seek(Pos)
84            x = []
85            y = []
86            w = []
87            S = File.readline()
88            while S and S[:4] != 'BANK' and S[0] != '#':
89                vals = S.split()
90                x.append(float(vals[0])/100.)               #CW: from centidegrees to degrees
91                f = float(vals[1])
92                if f > 0.0:
93                    y.append(f)
94                    w.append(1.0/f)
95                else:             
96                    y.append(0.0)
97                    w.append(0.0)
98                S = File.readline()
99            N = len(x)
100            return [np.array(x),np.array(y),np.array(w),np.zeros(N),np.zeros(N),np.zeros(N)]
101           
102        def GetESDdata(File,Pos,Bank):
103            File.seek(Pos)
104            cons = Bank.split()
105            if self.clockWd:
106                start = 0
107                step = 1
108            else:
109                start = float(cons[5])/100.0               #CW: from centidegrees to degrees
110                step = float(cons[6])/100.0
111            x = []
112            y = []
113            w = []
114            S = File.readline()
115            j = 0
116            while S and S[:4] != 'BANK' and S[0] != '#':
117                for i in range(0,80,16):
118                    xi = start+step*j
119                    yi = sfloat(S[i:i+8])
120                    ei = sfloat(S[i+8:i+16])
121                    x.append(xi)
122                    if yi > 0.0:
123                        y.append(yi)
124                        w.append(1.0/ei**2)
125                    else:             
126                        y.append(0.0)
127                        w.append(0.0)
128                    j += 1
129                S = File.readline()
130            N = len(x)
131            if self.clockWd:
132                x = Tmap2TOF(self.TimeMap,clockWd)
133            return [np.array(x),np.array(y),np.array(w),np.zeros(N),np.zeros(N),np.zeros(N)]
134       
135        def GetSTDdata(File,Pos,Bank):
136            File.seek(Pos)
137            cons = Bank.split()
138            Nch = int(cons[2])
139            if self.clockWd:
140                start = 0
141                step = 1
142            else:
143                start = float(cons[5])/100.0               #CW: from centidegrees to degrees
144                step = float(cons[6])/100.0                 #NB TOF 0.1*ms!
145            x = []
146            y = []
147            w = []
148            S = File.readline()
149            j = 0
150            while S and S[:4] != 'BANK' and S[0] != '#':
151                for i in range(0,80,8):
152                    xi = start+step*j
153                    ni = max(sint(S[i:i+2]),1)
154                    yi = max(sfloat(S[i+2:i+8]),0.0)
155                    if yi:
156                        vi = yi/ni
157                    else:
158                        yi = 0.0
159                        vi = 0.0
160                    j += 1
161                    if j < Nch:
162                        x.append(xi)
163                        if vi <= 0.:
164                            y.append(0.)
165                            w.append(0.)
166                        else:
167                            y.append(yi)
168                            w.append(1.0/vi)
169                S = File.readline()
170            N = len(x)
171            if self.clockWd:
172                x = Tmap2TOF(self.TimeMap,self.clockWd)[:-2]
173            return [np.array(x),np.array(y),np.array(w),np.zeros(N),np.zeros(N),np.zeros(N)]
174           
175        def GetALTdata(File,Pos,Bank):
176            File.seek(Pos)
177            cons = Bank.split()
178            x = []
179            y = []
180            w = []
181            S = File.readline()
182            j = 0
183            while S and S[:4] != 'BANK' and S[0] != '#':
184                for i in range(0,80,20):
185                    xi = sfloat(S[i:i+9])/3200.
186                    yi = sfloat(S[i+9:i+16])/1000.
187                    ei = sfloat(S[i+16:i+21])/1000.
188                    x.append(xi)
189                    if yi > 0.0:
190                        y.append(yi)
191                        w.append(1.0/ei**2)
192                    else:             
193                        y.append(0.0)
194                        w.append(0.0)
195                    j += 1
196                S = File.readline()
197            N = len(x)
198            if self.clockWd:
199                x = Tmap2TOF(self.TimeMap,clockWd)
200            return [np.array(x),np.array(y),np.array(w),np.zeros(N),np.zeros(N),np.zeros(N)]
201           
202        def GetTimeMap(File,Pos,TimeMap):
203            File.seek(Pos)
204            cons = TimeMap.split()
205            Nch = int(cons[1])
206            Nrec = int(cons[2])
207            clockWd = float(cons[4])/1000.          #in mus
208            TMap = np.zeros(Nch+2,dtype=int)
209            ind = 0
210            for i in range(Nrec):
211                S = File.readline().rstrip('\n')
212                vals = S.split()
213                for val in vals:
214                    TMap[ind] = int(val)
215                    ind += 1
216            TMap = np.reshape(TMap,(-1,3))
217            TMax = TMap[-1][0]
218            Nch = TMap[-2][0]+(TMax-TMap[-2][1]+TMap[-2][2]-1)/TMap[-2][2]
219            TMap[-1] = [Nch+1,TMax,0]
220            TMap = TMap.T
221            TMap[0] -= 1
222            return TMap.T,clockWd
223           
224        def Tmap2TOF(TMap,clockWd):
225            TOF = []
226            chWdt = []
227            Tch,T,Step = TMap[0]
228            for tmap in TMap[1:]:
229                tch,t,step = tmap
230                TOF += [T+Step*(i-Tch) for i in range(Tch,tch)]
231                Tch,T,Step = tmap
232            TOF = np.array(TOF)*clockWd
233            return TOF
234
235        x = []
236        y = []
237        w = []
238        Banks = []
239        Pos = []
240        rdbuffer = kwarg.get('buffer')
241        title = ''
242        comments = None
243#        selections = None
244
245        # reload previously saved values
246        if self.repeat and rdbuffer is not None:
247            Banks = rdbuffer.get('Banks')
248            Pos = rdbuffer.get('Pos')
249            self.selections = rdbuffer.get('selections')
250            comments = rdbuffer.get('comments')
251
252        # read through the file and find the beginning of each bank
253        # Save the offset (Pos), BANK line (Banks), comments for each bank
254        #
255        # This is going to need a fair amount of work to track line numbers
256        # in the input file.
257        if len(Banks) != len(Pos) or len(Banks) == 0:
258            try:
259                i = -1
260                while True:
261                    i += 1
262                    S = filepointer.readline()
263                    if len(S) == 0: break
264                       
265                    if i==0: # first line is always a comment
266                        self.errors = 'Error reading title'
267                        title = S[:-1]
268                        comments = [[title,]]
269                        continue
270                    if i==1 and S[:4].lower() == 'inst' and ':' in S:
271                        # 2nd line is instrument parameter file (optional)
272                        self.errors = 'Error reading instrument parameter filename'
273                        self.instparm = S.split(':')[1].strip('[]').strip()
274                        continue
275                    if S[0] == '#': # allow comments anywhere in the file
276                        # comments in fact should only preceed BANK lines
277                        comments[-1].append(S[:-1])
278                        continue
279                    if S[:4] == 'BANK':
280                        self.errors = 'Error reading bank:'
281                        self.errors += '  '+str(S)
282                        comments.append([title,])
283                        Banks.append(S)
284                        Pos.append(filepointer.tell())
285                    if S[:8] == 'TIME_MAP':
286                        if len(Banks) == 0:
287                            self.errors = 'Error reading time map before any bank lines'
288                        else:
289                            self.errors = 'Error reading time map after bank:\n  '+str(Banks[-1])
290                        self.TimeMap,self.clockWd = GetTimeMap(filepointer,filepointer.tell(),S)
291            except Exception as detail:
292                self.errors += '\n  '+str(detail)
293                print self.formatName+' scan error:'+str(detail) # for testing
294                import traceback
295                traceback.print_exc(file=sys.stdout)
296                return False
297
298        # Now select the bank to read
299        if not Banks: # use of ContentsValidator should prevent this error
300            print self.formatName+' scan error: no BANK records'
301            selblk = None # no block to choose
302            self.errors = 'No BANK records found (strange!)'
303            return False
304        elif len(Banks) == 1: # only one Bank, don't ask
305            selblk = 0
306        elif self.repeat and self.selections is not None:
307            # we were called to repeat the read
308            #print 'debug: repeat #',self.repeatcount,'selection',self.selections[self.repeatcount]
309            selblk = self.selections[self.repeatcount]
310            self.repeatcount += 1
311            if self.repeatcount >= len(self.selections): self.repeat = False
312        else:                       # choose from options
313            if not len(self.selections):    #use previous selection, otherwise...
314                self.selections = self.MultipleBlockSelector(
315                    Banks,
316                    ParentFrame=ParentFrame,
317                    title='Select Bank(s) to read from the list below',
318                    size=(600,100),
319                    header='Dataset Selector')
320            if len(self.selections) == 0: return False
321            selblk = self.selections[0] # select first in list
322            if len(self.selections) > 1: # prepare to loop through again
323                self.repeat = True
324                self.repeatcount = 1
325                if rdbuffer is not None:
326                    rdbuffer['Banks'] = Banks
327                    rdbuffer['Pos'] = Pos
328                    rdbuffer['selections'] = self.selections
329                    rdbuffer['comments'] = comments
330
331        # got a selection, now read it
332        Bank = Banks[selblk]
333        try:
334            bnkNo = int(Bank.split()[1])
335        except ValueError:
336            bnkNo = 1
337        try:
338            if 'FXYE' in Bank:
339                self.errors = 'Error reading FXYE data in Bank\n  '+Banks[selblk]
340                self.powderdata = GetFXYEdata(filepointer,Pos[selblk],Banks[selblk])
341            elif 'FXY' in Bank:
342                self.errors = 'Error reading FXY data in Bank\n  '+Banks[selblk]
343                self.powderdata = GetFXYdata(filepointer,Pos[selblk],Banks[selblk])
344            elif 'ESD' in Bank:
345                self.errors = 'Error reading ESD data in Bank\n  '+Banks[selblk]
346                self.powderdata = GetESDdata(filepointer,Pos[selblk],Banks[selblk])
347            elif 'STD' in Bank:
348                self.errors = 'Error reading STD data in Bank\n  '+Banks[selblk]
349                self.powderdata = GetSTDdata(filepointer,Pos[selblk],Banks[selblk])
350            elif 'ALT' in Bank:
351                self.errors = 'Error reading ALT data in Bank\n  '+Banks[selblk]
352                self.powderdata = GetALTdata(filepointer,Pos[selblk],Banks[selblk])
353            else:
354                self.errors = 'Error reading STD data in Bank\n  '+Banks[selblk]
355                self.powderdata = GetSTDdata(filepointer,Pos[selblk],Banks[selblk])
356        except Exception as detail:
357            self.errors += '\n  '+str(detail)
358            print self.formatName+' read error:'+str(detail) # for testing
359            import traceback
360            traceback.print_exc(file=sys.stdout)
361            return False
362
363        self.errors = 'Error processing information after read complete'
364        if comments is not None:
365            self.comments = comments[selblk]
366        self.powderentry[0] = filename
367        self.powderentry[1] = Pos # position offset (never used, I hope)
368        self.powderentry[2] = bnkNo #selblk+1 # bank number
369        self.idstring = ospath.basename(filename) + ' Bank '+str(bnkNo) #selblk+1)
370        self.numbanks=len(Banks)
371        # scan comments for temperature & radius
372        Temperature = 300.
373        for S in self.comments:
374            if 'Temp' in S.split('=')[0]:
375                try:
376                    Temperature = float(S.split('=')[1])
377                except:
378                    pass
379            elif 'Gonio' in S.split('=')[0]:
380                try:
381                    self.Sample['Gonio. radius'] = float(S.split('=')[1])
382                except:
383                    pass
384            elif 'Omega' in S.split('=')[0] or 'Theta' in S.split('=')[0]:  #HIPD weirdness
385                try:
386                    self.Sample['Omega'] = float(S.split('=')[1])
387                except:
388                    pass
389            elif 'Chi' in S.split('=')[0]:
390                try:
391                    self.Sample['Chi'] = float(S.split('=')[1])
392                except:
393                    pass                   
394            elif 'Phi' in S.split('=')[0]:
395                try:
396                    self.Sample['Phi'] = float(S.split('=')[1])
397                except:
398                    pass                   
399        self.Sample['Temperature'] = Temperature
400        return True       
401
402def sfloat(S):
403    'convert a string to a float, treating an all-blank string as zero'
404    if S.strip():
405        return float(S)
406    else:
407        return 0.0
408
409def sint(S):
410    'convert a string to an integer, treating an all-blank string as zero'
411    if S.strip():
412        return int(S)
413    else:
414        return 0
Note: See TracBrowser for help on using the repository browser.