source: trunk/GSASIIIO.py @ 246

Last change on this file since 246 was 246, checked in by vondreele, 11 years ago

fix copy controls for images - needed deepcopy
multi frame GE detector images now automatically summed when read.
also earlier commit - image calibration will compute a "suggested" new wavelength - not imposed however as it may be wrong.

  • Property svn:keywords set to Date Author Revision URL Id
File size: 46.1 KB
Line 
1"""GSASIIIO: functions for IO of data
2   Copyright: 2008, Robert B. Von Dreele (Argonne National Laboratory)
3"""
4########### SVN repository information ###################
5# $Date: 2011-02-14 00:20:15 +0000 (Mon, 14 Feb 2011) $
6# $Author: vondreele $
7# $Revision: 246 $
8# $URL: trunk/GSASIIIO.py $
9# $Id: GSASIIIO.py 246 2011-02-14 00:20:15Z vondreele $
10########### SVN repository information ###################
11import wx
12import math
13import numpy as np
14import cPickle
15import sys
16import random as ran
17import GSASIIpath
18import GSASIIgrid as G2gd
19import GSASIIspc as G2spc
20import GSASIIlattice as G2lat
21import GSASIIElem as G2el
22import os.path as ospath
23
24def sfloat(S):
25    if S.strip():
26        return float(S)
27    else:
28        return 0.0
29
30def sint(S):
31    if S.strip():
32        return int(S)
33    else:
34        return 0
35
36def SelectPowderData(self, filename):
37    """Selects banks of data from a filename of any GSAS powder data format
38    Input - filename: any GSAS powder data formatted file (currently STD, FXYE, FXY & ESD)
39    Returns - a list of banks to be read; each entry in list is a tuple containing:
40    filename: same as input filename
41    Pos: position for start of data; record just after BANK record
42    Bank: the BANK record
43    """
44    File = open(filename,'Ur')
45    Title = '''
46First line of this file:
47'''+File.readline()
48    dlg = wx.MessageDialog(self, Title, 'Is this the file you want?', 
49    wx.YES_NO | wx.ICON_QUESTION)
50    try:
51        result = dlg.ShowModal()
52    finally:
53        dlg.Destroy()
54    if result == wx.ID_NO: return (0,0)
55    Temperature = 300
56       
57    self.IparmName = GetInstrumentFile(self,filename)
58    if self.IparmName:
59        Iparm = GetInstrumentData(self.IparmName)
60    else:
61        Iparm = {}                                               #Assume CuKa lab data if no iparm file
62        Iparm['INS   HTYPE '] = 'PXC '
63        Iparm['INS  1 ICONS'] = '  1.540500  1.544300       0.0         0       0.7    0       0.5   '
64        Iparm['INS  1PRCF1 '] = '    3    8      0.01                                                '
65        Iparm['INS  1PRCF11'] = '   2.000000E+00  -2.000000E+00   5.000000E+00   0.000000E+00        '
66        Iparm['INS  1PRCF12'] = '   0.000000E+00   0.000000E+00   0.150000E-01   0.150000E-01        '
67    S = 1
68    Banks = []
69    Pos = []
70    FoundData = []
71    Comments = []
72    wx.BeginBusyCursor()
73    try:
74        while S:
75            S = File.readline()
76            if S[:1] != '#':
77                if S[:4] == 'BANK':
78                    Banks.append(S)
79                    Pos.append(File.tell())
80            else:
81                Comments.append(S[:-1])
82                if 'Temp' in S:
83                    Temperature = float(S[:-1].split()[-1])
84        File.close()
85    finally:
86        wx.EndBusyCursor()
87    if Comments:
88       print 'Comments on file:'
89       for Comment in Comments: print Comment
90    if Banks:
91        result = [0]
92        if len(Banks) >= 2:
93            dlg = wx.MultiChoiceDialog(self, 'Which scans do you want?', 'Select scans', Banks, wx.CHOICEDLG_STYLE)
94            try:
95                if dlg.ShowModal() == wx.ID_OK:
96                    result = dlg.GetSelections()
97                else:
98                    result = []
99            finally:
100                dlg.Destroy()
101        for i in result:
102            FoundData.append((filename,Pos[i],Banks[i]))
103    else:
104        dlg = wx.MessageDialog(self, 'ERROR - this is not a GSAS powder data file', 'No BANK records', wx.OK | wx.ICON_ERROR)
105        try:
106            result = dlg.ShowModal()
107        finally:
108            dlg.Destroy()
109    return FoundData,Iparm,Comments,Temperature
110
111def GetInstrumentFile(self,filename):
112    import os.path as op
113    dlg = wx.FileDialog(self,'Choose an instrument file','.', '', 'GSAS iparm file (*.prm)|*.prm|All files(*.*)|*.*', wx.OPEN)
114    if self.dirname: 
115        dlg.SetDirectory(self.dirname)
116        Tname = filename[:filename.index('.')]+'.prm'
117        if op.exists(Tname):
118            self.IparmName = Tname       
119    if self.IparmName: dlg.SetFilename(self.IparmName)
120    filename = ''
121    try:
122        if dlg.ShowModal() == wx.ID_OK:
123            filename = dlg.GetPath()
124    finally:
125        dlg.Destroy()
126    return filename
127
128def GetInstrumentData(IparmName):
129    file = open(IparmName, 'Ur')
130    S = 1
131    Iparm = {}
132    while S:
133        S = file.readline()
134        Iparm[S[:12]] = S[12:-1]
135    return Iparm
136   
137def GetPowderPeaks(fileName):
138    sind = lambda x: math.sin(x*math.pi/180.)
139    asind = lambda x: 180.*math.asin(x)/math.pi
140    Cuka = 1.54052
141    File = open(fileName,'Ur')
142    Comments = []
143    peaks = []
144    S = File.readline()
145    while S:
146        if S[:1] == '#':
147            Comments.append(S[:-1])
148        else:
149            item = S.split()
150            if len(item) == 1:
151                peaks.append([float(item[0]),1.0])
152            elif len(item) > 1:
153                peaks.append([float(item[0]),float(item[0])])
154        S = File.readline()
155    File.close()
156    if Comments:
157       print 'Comments on file:'
158       for Comment in Comments: print Comment
159    Peaks = []
160    if peaks[0][0] > peaks[-1][0]:          # d-spacings - assume CuKa
161        for peak in peaks:
162            dsp = peak[0]
163            sth = Cuka/(2.0*dsp)
164            if sth < 1.0:
165                tth = 2.0*asind(sth)
166            else:
167                break
168            Peaks.append([tth,peak[1],True,False,0,0,0,dsp,0.0])
169    else:                                   #2-thetas - assume Cuka (for now)
170        for peak in peaks:
171            tth = peak[0]
172            dsp = Cuka/(2.0*sind(tth/2.0))
173            Peaks.append([tth,peak[1],True,False,0,0,0,dsp,0.0])
174    return Comments,Peaks
175
176def GetPawleyPeaks(filename):
177    rt2ln2x2 = 2.35482
178    File = open(filename,'Ur')
179    PawleyPeaks = []
180    S = File.readline()         #skip header
181    S = File.readline()
182    item = S.split()
183    while S:
184        h,k,l = int(item[0]),int(item[1]),int(item[2])
185        mult = int(item[3])
186        tth = float(item[5])
187        sig = float(item[6])/rt2ln2x2
188        Iobs = float(item[7])*mult
189        PawleyPeaks.append([h,k,l,mult,tth,sig,False,Iobs,0.0,[]])
190        S = File.readline()
191        item = S.split()
192        if item[3] == '-100.0000':       #find trailer
193            break
194    File.close()
195    return PawleyPeaks
196   
197def GetHKLData(filename):
198    print 'Reading: '+filename
199    File = open(filename,'Ur')
200    HKLref = []
201    HKLmin = [1000,1000,1000]
202    HKLmax = [0,0,0]
203    FoMax = 0
204    ifFc = False
205    S = File.readline()
206    while '#' in S[0]:        #get past comments if any
207        S = File.readline()       
208    if '_' in S:         #cif style .hkl file
209        while 'loop_' not in S:         #skip preliminaries if any - can't have 'loop_' in them!
210            S = File.readline()       
211        S = File.readline()             #get past 'loop_' line
212        pos = 0
213        hpos = kpos = lpos = Fosqpos = Fcsqpos = sigpos = -1
214        while S:
215            if '_' in S:
216                if 'index_h' in S:
217                    hpos = pos
218                elif 'index_k' in S:
219                    kpos = pos
220                elif 'index_l' in S:
221                    lpos = pos
222                elif 'F_squared_meas' in S:
223                    Fosqpos = pos
224                elif 'F_squared_calc' in S:
225                    Fcsqpos = pos
226                elif 'F_squared_sigma' in S:
227                    sigpos = pos
228                pos += 1
229            else:
230                data = S.split()
231                if data:                    #avoid blank lines
232                    HKL = np.array([int(data[hpos]),int(data[kpos]),int(data[lpos])])
233                    h,k,l = HKL
234                    Fosq = float(data[Fosqpos])
235                    if sigpos != -1:
236                        sigFosq = float(data[sigpos])
237                    else:
238                        sigFosq = 1.
239                    if Fcsqpos != -1:
240                        Fcsq = float(data[Fcsqpos])
241                        if Fcsq:
242                            ifFc = True
243                    else:
244                        Fcsq = 0.
245                       
246                    HKLmin = [min(h,HKLmin[0]),min(k,HKLmin[1]),min(l,HKLmin[2])]
247                    HKLmax = [max(h,HKLmax[0]),max(k,HKLmax[1]),max(l,HKLmax[2])]
248                    FoMax = max(FoMax,Fosq)
249                    HKLref.append([HKL,Fosq,sigFosq,Fcsq,0,0,0])                 #room for Fcp, Fcpp & phase
250            S = File.readline()
251    else:                   #dumb h,k,l,Fo,sigFo .hkl file
252        while S:
253            h,k,l,Fo,sigFo = S.split()
254            HKL = np.array([int(h),int(k),int(l)])
255            h,k,l = HKL
256            Fo = float(Fo)
257            sigFo = float(sigFo)
258            HKLmin = [min(h,HKLmin[0]),min(k,HKLmin[1]),min(l,HKLmin[2])]
259            HKLmax = [max(h,HKLmax[0]),max(k,HKLmax[1]),max(l,HKLmax[2])]
260            FoMax = max(FoMax,Fo)
261            HKLref.append([HKL,Fo**2,2.*Fo*sigFo,0,0,0,0])                 #room for Fc, Fcp, Fcpp & phase
262            S = File.readline()
263    File.close()
264    return HKLref,HKLmin,HKLmax,FoMax,ifFc
265
266def GetPowderData(filename,Pos,Bank,DataType):
267    '''Reads one BANK of data from GSAS raw powder data file
268    input:
269    filename: GSAS raw powder file dataname
270    Pos: start of data in file just after BANK record
271    Bank: the BANK record
272    DataType: powder data type, e.g. "PXC" for Powder X-ray CW data
273    returns: list [x,y,e,yc,yb]
274    x: np.array of x-axis values
275    y: np.array of powder pattern intensities
276    w: np.array of w=sig(intensity)^2 values
277    yc: np.array of calc. intensities (zero)
278    yb: np.array of calc. background (zero)
279    yd: np.array of obs-calc profiles
280    '''
281    print 'Reading: '+filename
282    print 'Bank:    '+Bank[:-1]
283    if 'FXYE' in Bank:
284        return GetFXYEdata(filename,Pos,Bank,DataType)
285    elif 'FXY' in Bank:
286        return GetFXYdata(filename,Pos,Bank,DataType)
287    elif 'ESD' in Bank:
288        return GetESDdata(filename,Pos,Bank,DataType)
289    elif 'STD' in Bank:
290        return GetSTDdata(filename,Pos,Bank,DataType)
291    else:
292        return GetSTDdata(filename,Pos,Bank,DataType)
293    return []
294
295def GetFXYEdata(filename,Pos,Bank,DataType):
296    File = open(filename,'Ur')
297    File.seek(Pos)
298    x = []
299    y = []
300    w = []
301    S = File.readline()
302    while S and S[:4] != 'BANK':
303        vals = S.split()
304        if DataType[2] == 'C':
305            x.append(float(vals[0])/100.)               #CW: from centidegrees to degrees
306        elif DataType[2] == 'T':
307            x.append(float(vals[0])/1000.0)             #TOF: from musec to millisec
308        f = float(vals[1])
309        if f <= 0.0:
310            y.append(0.0)
311            w.append(1.0)
312        else:
313            y.append(float(vals[1]))
314            w.append(1.0/float(vals[2])**2)
315        S = File.readline()
316    File.close()
317    N = len(x)
318    return [np.array(x),np.array(y),np.array(w),np.zeros(N),np.zeros(N),np.zeros(N)]
319   
320def GetFXYdata(filename,Pos,Bank,DataType):
321    File = open(filename,'Ur')
322    File.seek(Pos)
323    x = []
324    y = []
325    w = []
326    S = File.readline()
327    while S and S[:4] != 'BANK':
328        vals = S.split()
329        if DataType[2] == 'C':
330            x.append(float(vals[0])/100.)               #CW: from centidegrees to degrees
331        elif DataType[2] == 'T':
332            x.append(float(vals[0])/1000.0)             #TOF: from musec to millisec
333        f = float(vals[1])
334        if f > 0.0:
335            y.append(f)
336            w.append(1.0/f)
337        else:             
338            y.append(0.0)
339            w.append(1.0)
340        S = File.readline()
341    File.close()
342    N = len(x)
343    return [np.array(x),np.array(y),np.array(w),np.zeros(N),np.zeros(N),np.zeros(N)]
344   
345def GetESDdata(filename,Pos,Bank,DataType):
346    File = open(filename,'Ur')
347    cons = Bank.split()
348    if DataType[2] == 'C':
349        start = float(cons[5])/100.0               #CW: from centidegrees to degrees
350        step = float(cons[6])/100.0
351    elif DataType[2] == 'T':
352        start = float(cons[5])/1000.0              #TOF: from musec to millisec
353        step = float(cons[6])/1000.0
354    File.seek(Pos)
355    x = []
356    y = []
357    w = []
358    S = File.readline()
359    j = 0
360    while S and S[:4] != 'BANK':
361        for i in range(0,80,16):
362            xi = start+step*j
363            yi = sfloat(S[i:i+8])
364            ei = sfloat(S[i+8:i+16])
365            x.append(xi)
366            if yi > 0.0:
367                y.append(yi)
368                w.append(1.0/ei**2)
369            else:             
370                y.append(0.0)
371                w.append(1.0)
372            j += 1
373        S = File.readline()
374    File.close()
375    N = len(x)
376    return [np.array(x),np.array(y),np.array(w),np.zeros(N),np.zeros(N),np.zeros(N)]
377
378def GetSTDdata(filename,Pos,Bank,DataType):
379    File = open(filename,'Ur')
380    cons = Bank.split()
381    Nch = cons[2]
382    if DataType[2] == 'C':
383        start = float(cons[5])/100.0               #CW: from centidegrees to degrees
384        step = float(cons[6])/100.0
385    elif DataType[2] == 'T':
386        start = float(cons[5])/1000.0              #TOF: from musec to millisec - not likely!
387        step = float(cons[6])/1000.0
388    File.seek(Pos)
389    x = []
390    y = []
391    w = []
392    S = File.readline()
393    j = 0
394    while S and S[:4] != 'BANK':
395        for i in range(0,80,8):
396            xi = start+step*j
397            ni = max(sint(S[i:i+2]),1)
398            yi = max(sfloat(S[i+2:i+8]),0.0)
399            if yi:
400                ei = math.sqrt(yi*ni)
401            else:
402                yi = 0.0
403                ei = 1.0
404            j += 1
405            if j < Nch:
406                x.append(xi)
407                y.append(yi)
408                w.append(1.0/ei**2)
409        S = File.readline()
410    File.close()
411    N = len(x)
412    return [np.array(x),np.array(y),np.array(w),np.zeros(N),np.zeros(N),np.zeros(N)]
413   
414def CheckImageFile(self,imagefile):
415    if not ospath.exists(imagefile):
416        dlg = wx.FileDialog(self, 'Bad image file name; choose name', '.', '',\
417        'Any image file (*.tif;*.tiff;*.mar*;*.avg;*.sum;*.img)\
418        |*.tif;*.tiff;*.mar*;*.avg;*.sum;*.img|\
419        Any detector tif (*.tif;*.tiff)|*.tif;*.tiff|\
420        MAR file (*.mar*)|*.mar*|\
421        GE Image (*.avg;*.sum)|*.avg;*.sum|\
422        ADSC Image (*.img)|*.img|\
423        All files (*.*)|*.*',wx.OPEN)
424        if self.dirname:
425            dlg.SetDirectory(self.dirname)
426        try:
427            dlg.SetFilename(ospath.split(imagefile)[1])
428            if dlg.ShowModal() == wx.ID_OK:
429                self.dirname = dlg.GetDirectory()
430                imagefile = dlg.GetPath()
431            else:
432                imagefile = False
433        finally:
434            dlg.Destroy()
435    return imagefile
436       
437def GetImageData(self,imagefile,imageOnly=False):       
438    ext = ospath.splitext(imagefile)[1]
439    Comments = []
440    if ext == '.tif':
441        Comments,Data,Npix,Image = GetTifData(imagefile)
442    elif ext == '.img':
443        Comments,Data,Npix,Image = GetImgData(imagefile)
444        Image[0][0] = 0
445    elif ext == '.mar3450' or ext == '.mar2300':
446        Comments,Data,Npix,Image = GetMAR345Data(imagefile)
447    elif ext in ['.sum','.avg','']:
448        Comments,Data,Npix,Image = GetGEsumData(imagefile)
449    elif ext == '.G2img':
450        return GetG2Image(imagefile)
451    if imageOnly:
452        return Image
453    else:
454        return Comments,Data,Npix,Image
455       
456def PutG2Image(filename,image):
457    File = open(filename,'wb')
458    cPickle.dump(image,File,1)
459    File.close()
460    return
461   
462def GetG2Image(filename):
463    File = open(filename,'rb')
464    image = cPickle.load(File)
465    File.close()
466    return image
467   
468def GetGEsumData(filename,imageOnly=False):
469    import struct as st
470    import array as ar
471    if not imageOnly:
472        print 'Read GE sum file: ',filename   
473    File = open(filename,'rb')
474    if '.sum' in filename:
475        head = ['GE detector sum data from APS 1-ID',]
476        sizexy = [2048,2048]
477    elif '.avg' in filename:
478        head = ['GE detector avg data from APS 1-ID',]
479        sizexy = [2048,2048]
480    else:
481        head = ['GE detector raw data from APS 1-ID',]
482        File.seek(18)
483        size,nframes = st.unpack('<ih',File.read(6))
484        sizexy = [2048,2048]
485        pos = 8192
486        File.seek(pos)
487    Npix = sizexy[0]*sizexy[1]
488    if '.sum' in filename:
489        image = np.array(ar.array('f',File.read(4*Npix)),dtype=np.int32)
490    elif '.avg' in filename:
491        image = np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)
492    else:
493        image = np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)
494        while nframes > 1:
495            image += np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)
496            nframes -= 1
497    image = np.reshape(image,(sizexy[1],sizexy[0]))
498    data = {'pixelSize':(200,200),'wavelength':0.15,'distance':250.0,'center':[204.8,204.8],'size':sizexy} 
499    File.close()   
500    if imageOnly:
501        return image
502    else:
503        return head,data,Npix,image
504       
505def GetImgData(filename,imageOnly=False):
506    import struct as st
507    import array as ar
508    if not imageOnly:
509        print 'Read ADSC img file: ',filename
510    File = open(filename,'rb')
511    head = File.read(511)
512    lines = head.split('\n')
513    head = []
514    center = [0,0]
515    for line in lines[1:-2]:
516        line = line.strip()[:-1]
517        if line:
518            if 'SIZE1' in line:
519                size = int(line.split('=')[1])
520                Npix = size*size
521            elif 'WAVELENGTH' in line:
522                wave = float(line.split('=')[1])
523            elif 'BIN' in line:
524                if line.split('=')[1] == '2x2':
525                    pixel=(102,102)
526                else:
527                    pixel = (51,51)
528            elif 'DISTANCE' in line:
529                distance = float(line.split('=')[1])
530            elif 'CENTER_X' in line:
531                center[0] = float(line.split('=')[1])
532            elif 'CENTER_Y' in line:
533                center[1] = float(line.split('=')[1])
534            head.append(line)
535    data = {'pixelSize':pixel,'wavelength':wave,'distance':distance,'center':center,'size':[size,size]}
536    image = []
537    row = 0
538    pos = 512
539    File.seek(pos)
540    image = np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)
541    image = np.reshape(image,(sizexy[1],sizexy[0]))
542#    image = np.zeros(shape=(size,size),dtype=np.int32)   
543#    while row < size:
544#        File.seek(pos)
545#        line = ar.array('H',File.read(2*size))
546#        image[row] = np.asarray(line)
547#        row += 1
548#        pos += 2*size
549    File.close()
550    if imageOnly:
551        return image
552    else:
553        return lines[1:-2],data,Npix,image
554       
555def GetMAR345Data(filename,imageOnly=False):
556    import array as ar
557    import struct as st
558    try:
559        import pack_f as pf
560    except:
561        msg = wx.MessageDialog(None, message="Unable to load the GSAS MAR image decompression, pack_f",
562                               caption="Import Error",
563                               style=wx.ICON_ERROR | wx.OK | wx.STAY_ON_TOP)
564        msg.ShowModal()
565        return None,None,None,None
566
567    if not imageOnly:
568        print 'Read Mar345 file: ',filename
569    File = open(filename,'rb')
570    head = File.read(4095)
571    numbers = st.unpack('<iiiiiiiiii',head[:40])
572    lines = head[128:].split('\n')
573    head = []
574    for line in lines:
575        line = line.strip()
576        if 'PIXEL' in line:
577            values = line.split()
578            pixel = (int(values[2]),int(values[4]))     #in microns
579        elif 'WAVELENGTH' in line:
580            wave = float(line.split()[1])
581        elif 'DISTANCE' in line:
582            distance = float(line.split()[1])           #in mm
583        elif 'CENTER' in line:
584            values = line.split()
585            center = [float(values[2])/10.,float(values[4])/10.]    #make in mm from pixels
586        if line: 
587            head.append(line)
588    data = {'pixelSize':pixel,'wavelength':wave,'distance':distance,'center':center}
589    for line in head:
590        if 'FORMAT' in line[0:6]:
591            items = line.split()
592            size = int(items[1])
593            Npix = size*size
594    pos = 4096
595    data['size'] = [size,size]
596    File.seek(pos)
597    line = File.read(8)
598    while 'CCP4' not in line:       #get past overflow list for now
599        line = File.read(8)
600        pos += 8
601    pos += 37
602    File.seek(pos)
603    raw = File.read()
604    File.close()
605    image = np.zeros(shape=(size,size),dtype=np.int32)
606    image = pf.pack_f(len(raw),raw,size,image)
607    if imageOnly:
608        return image.T              #transpose to get it right way around
609    else:
610        return head,data,Npix,image.T
611       
612def GetTifData(filename,imageOnly=False):
613    import struct as st
614    import array as ar
615    File = open(filename,'rb')
616    dataType = 5
617    try:
618        Meta = open(filename+'.metadata','Ur')
619        head = Meta.readlines()
620        for line in head:
621            line = line.strip()
622            if 'dataType' in line:
623                dataType = int(line.split('=')[1])
624        Meta.close()
625    except IOError:
626        print 'no metadata file found - will try to read file anyway'
627        head = ['no metadata file found',]
628       
629    tag = File.read(2)
630    byteOrd = '<'
631    if tag == 'II' and int(st.unpack('<h',File.read(2))[0]) == 42:     #little endian
632        IFD = int(st.unpack(byteOrd+'i',File.read(4))[0])
633    elif tag == 'MM' and int(st.unpack('>h',File.read(2))[0]) == 42:   #big endian
634        byteOrd = '>'
635        IFD = int(st.unpack(byteOrd+'i',File.read(4))[0])       
636    else:
637        lines = ['not a detector tiff file',]
638        return lines,0,0,0
639    File.seek(IFD)                                                  #get number of directory entries
640    NED = int(st.unpack(byteOrd+'h',File.read(2))[0])
641    IFD = {}
642    for ied in range(NED):
643        Tag,Type = st.unpack(byteOrd+'Hh',File.read(4))
644        nVal = st.unpack(byteOrd+'i',File.read(4))[0]
645        if Type == 1:
646            Value = st.unpack(byteOrd+nVal*'b',File.read(nVal))
647        elif Type == 2:
648            Value = st.unpack(byteOrd+'i',File.read(4))
649        elif Type == 3:
650            Value = st.unpack(byteOrd+nVal*'i',File.read(nVal*4))
651        elif Type == 4:
652            Value = st.unpack(byteOrd+nVal*'i',File.read(nVal*4))
653        elif Type == 5:
654            Value = st.unpack(byteOrd+nVal*'i',File.read(nVal*4))
655        elif Type == 11:
656            Value = st.unpack(byteOrd+nVal*'f',File.read(nVal*4))
657        IFD[Tag] = [Type,nVal,Value]
658    sizexy = [IFD[256][2][0],IFD[257][2][0]]
659    [nx,ny] = sizexy
660    Npix = nx*ny
661    if 272 in IFD:
662        ifd = IFD[272]
663        File.seek(ifd[2][0])
664        S = File.read(ifd[1])
665        if 'PILATUS' in S:
666            tifType = 'Pilatus'
667            dataType = 0
668            pixy = (172,172)
669            File.seek(4096)
670            if not imageOnly:
671                print 'Read Pilatus tiff file: ',filename
672            image = ar.array('L',File.read(4*Npix))
673            image = np.array(np.asarray(image),dtype=np.int32)
674    elif 262 in IFD and IFD[262][2][0] > 4:
675        tifType = 'DND'
676        pixy = (158,158)
677        File.seek(512)
678        if not imageOnly:
679            print 'Read DND SAX/WAX-detector tiff file: ',filename
680        image = np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)
681    elif sizexy == [1536,1536]:
682        tifType = 'APS Gold'
683        pixy = (150,150)
684        File.seek(64)
685        if not imageOnly:
686            print 'Read Gold tiff file:',filename
687        image = np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)
688    elif sizexy == [2048,2048] or sizexy == [1024,1024]:
689        if IFD[273][2][0] == 8:
690            if IFD[258][2][0] == 32:
691                tifType = 'PE'
692                pixy = (200,200)
693                File.seek(8)
694                if not imageOnly:
695                    print 'Read APS PE-detector tiff file: ',filename
696                if dataType == 5:
697                    image = np.array(ar.array('f',File.read(4*Npix)),dtype=np.int32)
698                else:
699                    image = np.array(ar.array('I',File.read(4*Npix)),dtype=np.int32)
700        elif IFD[273][2][0] == 4096:
701            tifType = 'MAR'
702            pixy = (158,158)
703            File.seek(4096)
704            if not imageOnly:
705                print 'Read MAR CCD tiff file: ',filename
706            image = np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)
707           
708    else:
709        lines = ['not a known detector tiff file',]
710        return lines,0,0,0
711       
712    image = np.reshape(image,(sizexy[1],sizexy[0]))
713    center = [pixy[0]*sizexy[0]/2000,pixy[1]*sizexy[1]/2000]
714    data = {'pixelSize':pixy,'wavelength':0.10,'distance':100.0,'center':center,'size':sizexy}
715    File.close()   
716    if imageOnly:
717        return image
718    else:
719        return head,data,Npix,image
720   
721def ProjFileOpen(self):
722    file = open(self.GSASprojectfile,'rb')
723    print 'load from file: ',self.GSASprojectfile
724    wx.BeginBusyCursor()
725    try:
726        while True:
727            try:
728                data = cPickle.load(file)
729            except EOFError:
730                break
731            datum = data[0]
732            print 'load: ',datum[0]
733           
734            Id = self.PatternTree.AppendItem(parent=self.root,text=datum[0])
735            if 'PWDR' in datum[0]:               
736                self.PatternTree.SetItemPyData(Id,datum[1][:3])     #temp. trim off junk
737            else:
738                self.PatternTree.SetItemPyData(Id,datum[1])
739            for datus in data[1:]:
740                print '    load: ',datus[0]
741                               
742                sub = self.PatternTree.AppendItem(Id,datus[0])
743                self.PatternTree.SetItemPyData(sub,datus[1])
744                                               
745            if 'IMG' in datum[0]:                   #retreive image default flag & data if set
746                Data = self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self,Id,'Image Controls'))
747                if Data['setDefault']:
748                    self.imageDefault = Data               
749        file.close()
750       
751    finally:
752        wx.EndBusyCursor()
753    print 'project load successful'
754    self.NewPlot = True
755   
756def ProjFileSave(self):
757    if not self.PatternTree.IsEmpty():
758        file = open(self.GSASprojectfile,'wb')
759        print 'save to file: ',self.GSASprojectfile
760        wx.BeginBusyCursor()
761        try:
762            item, cookie = self.PatternTree.GetFirstChild(self.root)
763            while item:
764                data = []
765                name = self.PatternTree.GetItemText(item)
766                print 'save: ',name
767                data.append([name,self.PatternTree.GetItemPyData(item)])
768                item2, cookie2 = self.PatternTree.GetFirstChild(item)
769                while item2:
770                    name = self.PatternTree.GetItemText(item2)
771                    print '    save: ',name
772                    data.append([name,self.PatternTree.GetItemPyData(item2)])
773                    item2, cookie2 = self.PatternTree.GetNextChild(item, cookie2)                           
774                item, cookie = self.PatternTree.GetNextChild(self.root, cookie)                           
775                cPickle.dump(data,file,1)
776            file.close()
777        finally:
778            wx.EndBusyCursor()
779        print 'project save successful'
780       
781def SaveIntegration(self,PickId,data):
782    azms = self.Integrate[1]
783    X = self.Integrate[2][:-1]
784    Xminmax = [X[0],X[-1]]
785    N = len(X)
786    Id = self.PatternTree.GetItemParent(PickId)
787    name = self.PatternTree.GetItemText(Id)
788    name = name.replace('IMG ','PWDR ')
789    Comments = self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self,Id, 'Comments'))
790    names = ['Type','Lam','Zero','Polariz.','U','V','W','X','Y','SH/L','Azimuth'] 
791    codes = [0 for i in range(11)]
792    parms = ['PXC',data['wavelength'],0.0,0.0,1.0,-1.0,0.3,0.0,1.0,0.0,0.0]
793    Azms = [(azms[i+1]+azms[i])/2. for i in range(len(azms)-1)]
794    for i,azm in enumerate(Azms):
795        item, cookie = self.PatternTree.GetFirstChild(self.root)
796        Id = 0
797        while item:
798            Name = self.PatternTree.GetItemText(item)
799            if name == Name:
800                Id = item
801            item, cookie = self.PatternTree.GetNextChild(self.root, cookie)
802        parms[10] = azm
803        Y = self.Integrate[0][i]
804        W = np.sqrt(Y)
805        Sample = {'Scale':[1.0,True],'Type':'Debye-Scherrer','Absorption':[0.0,False],'DisplaceX':[0.0,False],
806            'DisplaceY':[0.0,False],'Diffuse':[],'Temperature':300.,'Pressure':1.0,'Humidity':0.0,'Voltage':0.0,'Force':0.0}
807        if Id:
808            self.PatternTree.SetItemPyData(G2gd.GetPatternTreeItemId(self,Id, 'Comments'),Comments)                   
809            self.PatternTree.SetItemPyData(G2gd.GetPatternTreeItemId(self,Id,'Limits'),[tuple(Xminmax),Xminmax])
810            self.PatternTree.SetItemPyData(G2gd.GetPatternTreeItemId(self,Id,'Background'),[['chebyschev',1,3,1.0,0.0,0.0]])
811            self.PatternTree.SetItemPyData(G2gd.GetPatternTreeItemId(self,Id,'Instrument Parameters'),[tuple(parms),parms,codes,names])
812            self.PatternTree.SetItemPyData(G2gd.GetPatternTreeItemId(self,Id,'Peak List'),[])
813            self.PatternTree.SetItemPyData(G2gd.GetPatternTreeItemId(self,Id,'Index Peak List'),[])
814            self.PatternTree.SetItemPyData(G2gd.GetPatternTreeItemId(self,Id,'Unit Cells List'),[])             
815        else:
816            Id = self.PatternTree.AppendItem(parent=self.root,text=name+" Azm= %.2f"%(azm))
817            self.PatternTree.SetItemPyData(self.PatternTree.AppendItem(Id,text='Comments'),Comments)                   
818            self.PatternTree.SetItemPyData(self.PatternTree.AppendItem(Id,text='Limits'),[tuple(Xminmax),Xminmax])
819            self.PatternTree.SetItemPyData(self.PatternTree.AppendItem(Id,text='Background'),[['chebyschev',1,3,1.0,0.0,0.0]])
820            self.PatternTree.SetItemPyData(self.PatternTree.AppendItem(Id,text='Instrument Parameters'),[tuple(parms),parms,codes,names])
821            self.PatternTree.SetItemPyData(self.PatternTree.AppendItem(Id,text='Sample Parameters'),Sample)
822            self.PatternTree.SetItemPyData(self.PatternTree.AppendItem(Id,text='Peak List'),[])
823            self.PatternTree.SetItemPyData(self.PatternTree.AppendItem(Id,text='Index Peak List'),[])
824            self.PatternTree.SetItemPyData(self.PatternTree.AppendItem(Id,text='Unit Cells List'),[])             
825        self.PatternTree.SetItemPyData(Id,[[''],[np.array(X),np.array(Y),np.array(W),np.zeros(N),np.zeros(N),np.zeros(N)]])
826    self.PatternTree.SelectItem(Id)
827    self.PatternTree.Expand(Id)
828    self.PatternId = Id
829           
830def powderFxyeSave(self,powderfile):
831    file = open(powderfile,'w')
832    prm = open(powderfile.strip('fxye')+'prm','w')      #old style GSAS parm file
833    print 'save powder pattern to file: ',powderfile
834    wx.BeginBusyCursor()
835    Inst = self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self, \
836                    self.PickId, 'Instrument Parameters'))[1]
837    if len(Inst) == 11:             #single wavelength
838        lam1 = Inst[1]
839        lam2 = 0.0
840        GU,GV,GW = Inst[4:7]
841        LX,LY = Inst[7:9]
842        SL = HL = Inst[9]/2.0   
843    else:                           #Ka1 & Ka2
844        lam1 = Inst[1]
845        lam2 = Inst[2]
846        GU,GV,GW = Inst[6:9]
847        LX,LY = Inst[9:11]
848        SL = HL = Inst[11]/2.0   
849    prm.write( '            123456789012345678901234567890123456789012345678901234567890        '+'\n')
850    prm.write( 'INS   BANK      1                                                               '+'\n')
851    prm.write( 'INS   HTYPE   PXCR                                                              '+'\n')
852    prm.write(('INS  1 ICONS%10.7f%10.7f    0.0000               0.990    0     0.500   '+'\n')%(lam1,lam2))
853    prm.write( 'INS  1 IRAD     0                                                               '+'\n')
854    prm.write( 'INS  1I HEAD                                                                    '+'\n')
855    prm.write( 'INS  1I ITYP    0    0.0000  180.0000         1                                 '+'\n')
856    prm.write( 'INS  1PRCF1     3    8   0.00100                                                '+'\n')
857    prm.write(('INS  1PRCF11     %15.6g%15.6g%15.6g%15.6g   '+'\n')%(GU,GV,GW,0.0))
858    prm.write(('INS  1PRCF12     %15.6g%15.6g%15.6g%15.6g   '+'\n')%(LX,LY,SL,HL))
859    prm.close()
860    try:
861        x,y,w,yc,yb,yd = self.PatternTree.GetItemPyData(self.PickId)[1]
862        file.write(powderfile+'\n')
863        file.write('BANK 1 %d %d CONS %.2f %.2f 0 0 FXYE\n'%(len(x),len(x),\
864            100.*x[0],100.*(x[1]-x[0])))
865        s = list(np.sqrt(1./np.array(w)))       
866        XYW = zip(x,y,s)
867        for X,Y,S in XYW:
868            file.write("%15.6g %15.6g %15.6g\n" % (100.*X,Y,S))
869        file.close()
870    finally:
871        wx.EndBusyCursor()
872    print 'powder pattern file written'
873       
874def powderXyeSave(self,powderfile):
875    file = open(powderfile,'w')
876    print 'save powder pattern to file: ',powderfile
877    wx.BeginBusyCursor()
878    try:
879        x,y,w,yc,yb,yd = self.PatternTree.GetItemPyData(self.PickId)[1]
880        XYW = zip(x,y,w)
881        for X,Y,W in XYW:
882            file.write("%15.6g %15.6g %15.6g\n" % (X,Y,W))
883        file.close()
884    finally:
885        wx.EndBusyCursor()
886    print 'powder pattern file written'
887   
888def PeakListSave(self,file,peaks):
889    print 'save peak list to file: ',self.peaklistfile
890    if not peaks:
891        dlg = wx.MessageDialog(self, 'No peaks!', 'Nothing to save!', wx.OK)
892        try:
893            result = dlg.ShowModal()
894        finally:
895            dlg.Destroy()
896        return
897    for peak in peaks:
898        file.write("%10.4f %12.2f %10.3f %10.3f \n" % \
899            (peak[0],peak[2],peak[4],peak[6]))
900    print 'peak list saved'
901             
902def IndexPeakListSave(self,peaks):
903    file = open(self.peaklistfile,'wa')
904    print 'save index peak list to file: ',self.peaklistfile
905    wx.BeginBusyCursor()
906    try:
907        if not peaks:
908            dlg = wx.MessageDialog(self, 'No peaks!', 'Nothing to save!', wx.OK)
909            try:
910                result = dlg.ShowModal()
911            finally:
912                dlg.Destroy()
913            return
914        for peak in peaks:
915            file.write("%12.6f\n" % (peak[7]))
916        file.close()
917    finally:
918        wx.EndBusyCursor()
919    print 'index peak list saved'
920   
921def ReadEXPPhase(self,filename):
922    file = open(filename, 'Ur')
923    Phase = {}
924    S = 1
925    Expr = [{},{},{},{},{},{},{},{},{}]
926    while S:
927        S = file.readline()
928        if 'EXPR NPHAS' in S[:12]:
929            Num = S[12:-1].count('0')
930            NPhas = S[12:-1].split()
931        if 'CRS' in S[:3]:
932            N = int(S[3:4])-1
933            Expr[N][S[:12]] = S[12:-1]
934    file.close()
935    PNames = []
936    for n,N in enumerate(NPhas):
937        if N != '0':
938            result = n
939            key = 'CRS'+str(n+1)+'    PNAM'
940            PNames.append(Expr[n][key])
941    if Num < 8:
942        dlg = wx.SingleChoiceDialog(self, 'Which phase to read?', 'Read phase data', PNames, wx.CHOICEDLG_STYLE)
943        try:
944            if dlg.ShowModal() == wx.ID_OK:
945                result = dlg.GetSelection()
946        finally:
947            dlg.Destroy()       
948    EXPphase = Expr[result]
949    keyList = EXPphase.keys()
950    keyList.sort()
951    SGData = {}
952    if NPhas[result] == '1':
953        Ptype = 'nuclear'
954    elif NPhas[result] in ['2','3']:
955        Ptype = 'magnetic'
956    elif NPhas[result] == '4':
957        Ptype = 'macromolecular'
958    elif NPhas[result] == '10':
959        Ptype = 'Pawley'
960    for key in keyList:
961        if 'PNAM' in key:
962           PhaseName = EXPphase[key].strip()
963        elif 'ABC   ' in key:
964            abc = [float(EXPphase[key][:10]),float(EXPphase[key][10:20]),float(EXPphase[key][20:30])]                       
965        elif 'ANGLES' in key:
966            angles = [float(EXPphase[key][:10]),float(EXPphase[key][10:20]),float(EXPphase[key][20:30])]                                               
967        elif 'SG SYM' in key:
968            SpGrp = EXPphase[key][:15].strip()
969            E,SGData = G2spc.SpcGroup(SpGrp)
970    Atoms = []
971    if Ptype == 'nuclear':
972        for key in keyList:
973            if 'AT' in key:
974                if key[11:] == 'A':
975                    S = EXPphase[key]
976                elif key[11:] == 'B':
977                    S += EXPphase[key]
978                    Atom = [S[50:58].strip(),S[:10].strip(),'',
979                        float(S[10:20]),float(S[20:30]),float(S[30:40]),
980                        float(S[40:50]),'',int(S[60:62]),S[130:131]]
981                    if Atom[9] == 'I':
982                        Atom += [float(S[68:78]),0.,0.,0.,0.,0.,0.]
983                    elif Atom[9] == 'A':
984                        Atom += [0.0,float(S[68:78]),float(S[78:88]),
985                            float(S[88:98]),float(S[98:108]),
986                            float(S[108:118]),float(S[118:128])]
987                    XYZ = Atom[3:6]
988                    Atom[7],Atom[8] = G2spc.SytSym(XYZ,SGData)
989                    Atom.append(ran.randint(0,sys.maxint))
990                    Atoms.append(Atom)
991    elif Ptype == 'macromolecular':
992        for key in keyList:
993            if 'AT' in key[6:8]:
994                S = EXPphase[key]
995                Atom = [int(S[56:60]),S[50:54].strip().upper(),S[54:56],
996                    S[46:51].strip(),S[:8].strip(),'',
997                    float(S[16:24]),float(S[24:32]),float(S[32:40]),
998                    float(S[8:16]),'1',1,'I',float(S[40:46]),0,0,0,0,0,0]
999                XYZ = Atom[6:9]
1000                Atom[10],Atom[11] = G2spc.SytSym(XYZ,SGData)
1001                Atom.append(ran.randint(0,sys.maxint))
1002                Atoms.append(Atom)
1003    Volume = G2lat.calc_V(G2lat.cell2A(abc+angles))
1004    Phase['General'] = {'Name':PhaseName,'Type':Ptype,'SGData':SGData,
1005        'Cell':[False,]+abc+angles+[Volume,]}
1006    Phase['Atoms'] = Atoms
1007    Phase['Drawing'] = {}
1008    Phase['Histograms'] = {}
1009    return Phase
1010       
1011def ReadPDBPhase(filename):
1012    EightPiSq = 8.*math.pi**2
1013    file = open(filename, 'Ur')
1014    Phase = {}
1015    Title = ''
1016    Compnd = ''
1017    Atoms = []
1018    A = np.zeros(shape=(3,3))
1019    S = file.readline()
1020    while S:
1021        Atom = []
1022        if 'TITLE' in S[:5]:
1023            Title = S[10:72].strip()
1024            S = file.readline()
1025        elif 'COMPND    ' in S[:10]:
1026            Compnd = S[10:72].strip()
1027            S = file.readline()
1028        elif 'CRYST' in S[:5]:
1029            abc = S[7:34].split()
1030            angles = S[34:55].split()
1031            cell=[float(abc[0]),float(abc[1]),float(abc[2]),
1032                float(angles[0]),float(angles[1]),float(angles[2])]
1033            Volume = G2lat.calc_V(G2lat.cell2A(cell))
1034            AA,AB = G2lat.cell2AB(cell)
1035            SpGrp = S[55:65]
1036            E,SGData = G2spc.SpcGroup(SpGrp)
1037            if E: 
1038                print ' ERROR in space group symbol ',SpGrp,' in file ',filename
1039                print ' N.B.: make sure spaces separate axial fields in symbol' 
1040                print G2spc.SGErrors(E)
1041                return None
1042            SGlines = G2spc.SGPrint(SGData)
1043            for line in SGlines: print line
1044            S = file.readline()
1045        elif 'SCALE' in S[:5]:
1046            V = (S[10:41].split())
1047            A[int(S[5])-1] = [float(V[0]),float(V[1]),float(V[2])]
1048            S = file.readline()
1049        elif 'ATOM' in S[:4] or 'HETATM' in S[:6]:
1050            XYZ = [float(S[31:39]),float(S[39:47]),float(S[47:55])]
1051            XYZ = np.inner(AB,XYZ)
1052            SytSym,Mult = G2spc.SytSym(XYZ,SGData)
1053            Uiso = float(S[61:67])/EightPiSq
1054            Type = S[12:14].upper()
1055            if Type[0] in '123456789':
1056                Type = Type[1:]
1057            Atom = [S[22:27].strip(),S[17:20].upper(),S[20:22],
1058                S[12:17].strip(),Type.strip(),'',XYZ[0],XYZ[1],XYZ[2],
1059                float(S[55:61]),SytSym,Mult,'I',Uiso,0,0,0,0,0,0]
1060            S = file.readline()
1061            if 'ANISOU' in S[:6]:
1062                Uij = S[30:72].split()
1063                Uij = [float(Uij[0])/10000.,float(Uij[1])/10000.,float(Uij[2])/10000.,
1064                    float(Uij[3])/10000.,float(Uij[4])/10000.,float(Uij[5])/10000.]
1065                Atom = Atom[:14]+Uij
1066                Atom[12] = 'A'
1067                S = file.readline()
1068            Atom.append(ran.randint(0,sys.maxint))
1069            Atoms.append(Atom)
1070        else:           
1071            S = file.readline()
1072    file.close()
1073    if Title:
1074        PhaseName = Title
1075    elif Compnd:
1076        PhaseName = Compnd
1077    else:
1078        PhaseName = 'None'
1079    Phase['General'] = {'Name':PhaseName,'Type':'macromolecular','SGData':SGData,
1080        'Cell':[False,]+cell+[Volume,]}
1081    Phase['Atoms'] = Atoms
1082    Phase['Drawing'] = {}
1083    Phase['Histograms'] = {}
1084   
1085    return Phase
1086   
1087def ReadCIFPhase(filename):
1088    anisoNames = ['aniso_u_11','aniso_u_22','aniso_u_33','aniso_u_12','aniso_u_13','aniso_u_23']
1089    file = open(filename, 'Ur')
1090    Phase = {}
1091    Title = ospath.split(filename)[-1]
1092    print '\n Reading cif file: ',Title
1093    Compnd = ''
1094    Atoms = []
1095    A = np.zeros(shape=(3,3))
1096    S = file.readline()
1097    while S:
1098        if '_symmetry_space_group_name_H-M' in S:
1099            SpGrp = S.split("_symmetry_space_group_name_H-M")[1].strip().strip('"').strip("'")
1100            E,SGData = G2spc.SpcGroup(SpGrp)
1101            if E:
1102                print ' ERROR in space group symbol ',SpGrp,' in file ',filename
1103                print ' N.B.: make sure spaces separate axial fields in symbol' 
1104                print G2spc.SGErrors(E)
1105                return None
1106            S = file.readline()
1107        elif '_cell' in S:
1108            if '_cell_length_a' in S:
1109                a = S.split('_cell_length_a')[1].strip().strip('"').strip("'").split('(')[0]
1110            elif '_cell_length_b' in S:
1111                b = S.split('_cell_length_b')[1].strip().strip('"').strip("'").split('(')[0]
1112            elif '_cell_length_c' in S:
1113                c = S.split('_cell_length_c')[1].strip().strip('"').strip("'").split('(')[0]
1114            elif '_cell_angle_alpha' in S:
1115                alp = S.split('_cell_angle_alpha')[1].strip().strip('"').strip("'").split('(')[0]
1116            elif '_cell_angle_beta' in S:
1117                bet = S.split('_cell_angle_beta')[1].strip().strip('"').strip("'").split('(')[0]
1118            elif '_cell_angle_gamma' in S:
1119                gam = S.split('_cell_angle_gamma')[1].strip().strip('"').strip("'").split('(')[0]
1120            S = file.readline()
1121        elif 'loop_' in S:
1122            labels = {}
1123            i = 0
1124            while S:
1125                S = file.readline()
1126                if '_atom_site' in S.strip()[:10]:
1127                    labels[S.strip().split('_atom_site_')[1].lower()] = i
1128                    i += 1
1129                else:
1130                    break
1131            if labels:
1132                if 'aniso_label' not in labels:
1133                    while S:
1134                        atom = ['','','',0,0,0,0,'','','I',0.01,0,0,0,0,0,0]
1135                        S.strip()
1136                        if len(S.split()) != len(labels):
1137                            if 'loop_' in S:
1138                                break
1139                            S += file.readline().strip()
1140                        data = S.split()
1141                        if len(data) != len(labels):
1142                            break
1143                        for key in labels:
1144                            if key == 'type_symbol':
1145                                atom[1] = data[labels[key]]
1146                            elif key == 'label':
1147                                atom[0] = data[labels[key]]
1148                            elif key == 'fract_x':
1149                                atom[3] = float(data[labels[key]].split('(')[0])
1150                            elif key == 'fract_y':
1151                                atom[4] = float(data[labels[key]].split('(')[0])
1152                            elif key == 'fract_z':
1153                                atom[5] = float(data[labels[key]].split('(')[0])
1154                            elif key == 'occupancy':
1155                                atom[6] = float(data[labels[key]].split('(')[0])
1156                            elif key == 'thermal_displace_type':
1157                                if data[labels[key]].lower() == 'uiso':
1158                                    atom[9] = 'I'
1159                                    atom[10] = float(data[labels['u_iso_or_equiv']].split('(')[0])
1160                                else:
1161                                    atom[9] = 'A'
1162                                    atom[10] = 0.0
1163                                   
1164                        atom[7],atom[8] = G2spc.SytSym(atom[3:6],SGData)
1165                        atom.append(ran.randint(0,sys.maxint))
1166                        Atoms.append(atom)
1167                        S = file.readline()
1168                else:
1169                    while S:
1170                        S.strip()
1171                        data = S.split()
1172                        if len(data) != len(labels):
1173                            break
1174                        name = data[labels['aniso_label']]
1175                        for atom in Atoms:
1176                            if name == atom[0]:
1177                                for i,uname in enumerate(anisoNames):
1178                                    atom[i+11] = float(data[labels[uname]].split('(')[0])
1179                        S = file.readline()
1180                                                                       
1181        else:           
1182            S = file.readline()
1183    file.close()
1184    if Title:
1185        PhaseName = Title
1186    else:
1187        PhaseName = 'None'
1188    SGlines = G2spc.SGPrint(SGData)
1189    for line in SGlines: print line
1190    cell = [float(a),float(b),float(c),float(alp),float(bet),float(gam)]
1191    Volume = G2lat.calc_V(G2lat.cell2A(cell))
1192    Phase['General'] = {'Name':PhaseName,'Type':'nuclear','SGData':SGData,
1193        'Cell':[False,]+cell+[Volume,]}
1194    Phase['Atoms'] = Atoms
1195    Phase['Drawing'] = {}
1196    Phase['Histograms'] = {}
1197   
1198    return Phase
Note: See TracBrowser for help on using the repository browser.