source: trunk/imports/G2pwd_fxye.py @ 1515

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

allow '#' comments in GSAS iparm files.
add class G2ColumnIDDialog for loading sample data for all histograms from a text file. This file has optional '#' comments followed by columns of data items each line begins with a file name that matches the file names for the histograms.
This class allows simple arithmetic modification of a column of data.
G2MultiChoiceDialog now allows selection of a block of items
Add 'Time' to Sample parameters - typically clock time
Further work on modulated structures
change default .gsa GSAS powder data extension to .gda - a lot more common.

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