source: trunk/imports/G2pwd_fxye.py @ 2824

Last change on this file since 2824 was 2824, checked in by vondreele, 5 years ago

fix 2D strain d-zero calc
set Refine dlg=None default for easier scripting use
fix to put (?) space in BANK n

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