source: trunk/imports/G2pwd_fxye.py @ 2516

Last change on this file since 2516 was 2516, checked in by toby, 5 years ago

revise import to not assume Bank 1 with multibank instparm files; deal with unicode problem in CIF files; improve atoms use of selection from menu; add Pawley menu variable selection (plenty more to do)

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