source: trunk/imports/G2pwd_fxye.py @ 1767

Last change on this file since 1767 was 1767, checked in by vondreele, 8 years ago

fix texture calculation - now matches GSAS
  • 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: 2015-03-30 12:39:25 +0000 (Mon, 30 Mar 2015) $
4# $Author: vondreele $
5# $Revision: 1767 $
6# $URL: trunk/imports/G2pwd_fxye.py $
7# $Id: G2pwd_fxye.py 1767 2015-03-30 12:39:25Z 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: 1767 $")
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[8:].split()
205            mapNo = int(cons[0])
206            Nch = int(cons[1])
207            Nrec = int(cons[2])
208            clockWd = float(cons[4])/1000.          #in mus
209            TMap = np.zeros(Nch+2,dtype=int)
210            ind = 0
211            for i in range(Nrec):
212                S = File.readline().rstrip('\n')
213                vals = S.split()
214                for val in vals:
215                    TMap[ind] = int(val)
216                    ind += 1
217            TMap = np.reshape(TMap,(-1,3))
218            TMax = TMap[-1][0]
219            Nch = TMap[-2][0]+(TMax-TMap[-2][1]+TMap[-2][2]-1)/TMap[-2][2]
220            TMap[-1] = [Nch+1,TMax,0]
221            TMap = TMap.T
222            TMap[0] -= 1
223            return TMap.T,clockWd
224           
225        def Tmap2TOF(TMap,clockWd):
226            TOF = []
227            chWdt = []
228            Tch,T,Step = TMap[0]
229            for tmap in TMap[1:]:
230                tch,t,step = tmap
231                TOF += [T+Step*(i-Tch) for i in range(Tch,tch)]
232                Tch,T,Step = tmap
233            TOF = np.array(TOF)*clockWd
234            return TOF
235
236        x = []
237        y = []
238        w = []
239        Banks = []
240        Pos = []
241        rdbuffer = kwarg.get('buffer')
242        title = ''
243        comments = None
244#        selections = None
245
246        # reload previously saved values
247        if self.repeat and rdbuffer is not None:
248            Banks = rdbuffer.get('Banks')
249            Pos = rdbuffer.get('Pos')
250            self.selections = rdbuffer.get('selections')
251            comments = rdbuffer.get('comments')
252
253        # read through the file and find the beginning of each bank
254        # Save the offset (Pos), BANK line (Banks), comments for each bank
255        #
256        # This is going to need a fair amount of work to track line numbers
257        # in the input file.
258        if len(Banks) != len(Pos) or len(Banks) == 0:
259            try:
260                i = -1
261                while True:
262                    i += 1
263                    S = filepointer.readline()
264                    if len(S) == 0: break
265                       
266                    if i==0: # first line is always a comment
267                        self.errors = 'Error reading title'
268                        title = S[:-1]
269                        comments = [[title,]]
270                        continue
271                    if i==1 and S[:4].lower() == 'inst' and ':' in S:
272                        # 2nd line is instrument parameter file (optional)
273                        self.errors = 'Error reading instrument parameter filename'
274                        self.instparm = S.split(':')[1].strip('[]').strip()
275                        continue
276                    if S[0] == '#': # allow comments anywhere in the file
277                        # comments in fact should only preceed BANK lines
278                        comments[-1].append(S[:-1])
279                        continue
280                    if S[:4] == 'BANK':
281                        self.errors = 'Error reading bank:'
282                        self.errors += '  '+str(S)
283                        comments.append([title,])
284                        Banks.append(S)
285                        Pos.append(filepointer.tell())
286                    if S[:8] == 'TIME_MAP':     #assumes one time map; HIPPO has multiple time maps
287                        if len(Banks) == 0:
288                            self.errors = 'Error reading time map before any bank lines'
289                        else:
290                            self.errors = 'Error reading time map after bank:\n  '+str(Banks[-1])
291                        self.TimeMap,self.clockWd = GetTimeMap(filepointer,filepointer.tell(),S)
292            except Exception as detail:
293                self.errors += '\n  '+str(detail)
294                print self.formatName+' scan error:'+str(detail) # for testing
295                import traceback
296                traceback.print_exc(file=sys.stdout)
297                return False
298
299        # Now select the bank to read
300        if not Banks: # use of ContentsValidator should prevent this error
301            print self.formatName+' scan error: no BANK records'
302            selblk = None # no block to choose
303            self.errors = 'No BANK records found (strange!)'
304            return False
305        elif len(Banks) == 1: # only one Bank, don't ask
306            selblk = 0
307        elif self.repeat and self.selections is not None:
308            # we were called to repeat the read
309            #print 'debug: repeat #',self.repeatcount,'selection',self.selections[self.repeatcount]
310            selblk = self.selections[self.repeatcount]
311            self.repeatcount += 1
312            if self.repeatcount >= len(self.selections): self.repeat = False
313        else:                       # choose from options
314            if not len(self.selections):    #use previous selection, otherwise...
315                self.selections = self.MultipleBlockSelector(
316                    Banks,
317                    ParentFrame=ParentFrame,
318                    title='Select Bank(s) to read from the list below',
319                    size=(600,100),
320                    header='Dataset Selector')
321            if len(self.selections) == 0: return False
322            selblk = self.selections[0] # select first in list
323            if len(self.selections) > 1: # prepare to loop through again
324                self.repeat = True
325                self.repeatcount = 1
326                if rdbuffer is not None:
327                    rdbuffer['Banks'] = Banks
328                    rdbuffer['Pos'] = Pos
329                    rdbuffer['selections'] = self.selections
330                    rdbuffer['comments'] = comments
331
332        # got a selection, now read it
333        Bank = Banks[selblk]
334        try:
335            bnkNo = int(Bank.split()[1])
336        except ValueError:
337            bnkNo = 1
338        try:
339            if 'FXYE' in Bank:
340                self.errors = 'Error reading FXYE data in Bank\n  '+Banks[selblk]
341                self.powderdata = GetFXYEdata(filepointer,Pos[selblk],Banks[selblk])
342            elif 'FXY' in Bank:
343                self.errors = 'Error reading FXY data in Bank\n  '+Banks[selblk]
344                self.powderdata = GetFXYdata(filepointer,Pos[selblk],Banks[selblk])
345            elif 'ESD' in Bank:
346                self.errors = 'Error reading ESD data in Bank\n  '+Banks[selblk]
347                self.powderdata = GetESDdata(filepointer,Pos[selblk],Banks[selblk])
348            elif 'STD' in Bank:
349                self.errors = 'Error reading STD data in Bank\n  '+Banks[selblk]
350                self.powderdata = GetSTDdata(filepointer,Pos[selblk],Banks[selblk])
351            elif 'ALT' in Bank:
352                self.errors = 'Error reading ALT data in Bank\n  '+Banks[selblk]
353                self.powderdata = GetALTdata(filepointer,Pos[selblk],Banks[selblk])
354            else:
355                self.errors = 'Error reading STD data in Bank\n  '+Banks[selblk]
356                self.powderdata = GetSTDdata(filepointer,Pos[selblk],Banks[selblk])
357        except Exception as detail:
358            self.errors += '\n  '+str(detail)
359            print self.formatName+' read error:'+str(detail) # for testing
360            import traceback
361            traceback.print_exc(file=sys.stdout)
362            return False
363
364        self.errors = 'Error processing information after read complete'
365        if comments is not None:
366            self.comments = comments[selblk]
367        self.powderentry[0] = filename
368        self.powderentry[1] = Pos # position offset (never used, I hope)
369        self.powderentry[2] = bnkNo #selblk+1 # bank number
370        self.idstring = ospath.basename(filename) + ' Bank '+str(bnkNo) #selblk+1)
371        self.numbanks=len(Banks)
372        # scan comments for temperature & radius
373        Temperature = 300.
374        for S in self.comments:
375            if 'Temp' in S.split('=')[0]:
376                try:
377                    Temperature = float(S.split('=')[1])
378                except:
379                    pass
380            elif 'Gonio' in S.split('=')[0]:
381                try:
382                    self.Sample['Gonio. radius'] = float(S.split('=')[1])
383                except:
384                    pass
385            elif 'Omega' in S.split('=')[0] or 'Theta' in S.split('=')[0]:  #HIPD weirdness
386                try:
387                    self.Sample['Omega'] = 90.-float(S.split('=')[1])
388                except:
389                    pass
390            elif 'Chi' in S.split('=')[0]:
391                try:
392                    self.Sample['Chi'] = -float(S.split('=')[1])
393                except:
394                    pass                   
395            elif 'Phi' in S.split('=')[0]:
396                try:
397                    self.Sample['Phi'] = float(S.split('=')[1])
398                except:
399                    pass                   
400        self.Sample['Temperature'] = Temperature
401        return True       
402
403def sfloat(S):
404    'convert a string to a float, treating an all-blank string as zero'
405    if S.strip():
406        return float(S)
407    else:
408        return 0.0
409
410def sint(S):
411    'convert a string to an integer, treating an all-blank string as zero'
412    if S.strip():
413        return int(S)
414    else:
415        return 0
Note: See TracBrowser for help on using the repository browser.