source: trunk/imports/G2pwd_fxye.py @ 2152

Last change on this file since 2152 was 2152, checked in by vondreele, 6 years ago

all PWDR exporters will make file name from histogram name
allow read of multibank data
alert user to duplicate histograms (by name)
rename data will not change Bank or Azm part of histogram name
fix L&R plotting commands for TOF data

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