source: trunk/GSASIIIO.py @ 1138

Last change on this file since 1138 was 1138, checked in by toby, 8 years ago

major constraints revision

  • Property svn:eol-style set to native
  • Property svn:keywords set to Date Author Revision URL Id
File size: 93.5 KB
Line 
1# -*- coding: utf-8 -*-
2########### SVN repository information ###################
3# $Date: 2013-11-07 18:12:55 +0000 (Thu, 07 Nov 2013) $
4# $Author: toby $
5# $Revision: 1138 $
6# $URL: trunk/GSASIIIO.py $
7# $Id: GSASIIIO.py 1138 2013-11-07 18:12:55Z toby $
8########### SVN repository information ###################
9'''
10*GSASIIIO: Misc I/O routines*
11=============================
12
13Module with miscellaneous routines for input and output. Many
14are GUI routines to interact with user.
15
16Includes support for image reading.
17
18Also includes base classes for data import routines.
19
20'''
21"""GSASIIIO: functions for IO of data
22   Copyright: 2008, Robert B. Von Dreele (Argonne National Laboratory)
23"""
24import wx
25import math
26import numpy as np
27import cPickle
28import sys
29import re
30import random as ran
31import GSASIIpath
32GSASIIpath.SetVersionNumber("$Revision: 1138 $")
33import GSASIIgrid as G2gd
34import GSASIIspc as G2spc
35import GSASIIlattice as G2lat
36import GSASIIpwdGUI as G2pdG
37import GSASIIElem as G2el
38import GSASIIstrIO as G2stIO
39import GSASIImapvars as G2mv
40import os
41import os.path as ospath
42
43def sfloat(S):
44    'Convert a string to float. An empty field is treated as zero'
45    if S.strip():
46        return float(S)
47    else:
48        return 0.0
49
50def sint(S):
51    'Convert a string to int. An empty field is treated as zero'
52    if S.strip():
53        return int(S)
54    else:
55        return 0
56
57def trim(val):
58    '''Simplify a string containing leading and trailing spaces
59    as well as newlines, tabs, repeated spaces etc. into a shorter and
60    more simple string, by replacing all ranges of whitespace
61    characters with a single space.
62
63    :param str val: the string to be simplified
64
65    :returns: the (usually) shortened version of the string
66    '''
67    return re.sub('\s+', ' ', val).strip()
68
69def makeInstDict(names,data,codes):
70    inst = dict(zip(names,zip(data,data,codes)))
71    for item in inst:
72        inst[item] = list(inst[item])
73    return inst
74
75def FileDlgFixExt(dlg,file):
76    'this is needed to fix a problem in linux wx.FileDialog'
77    ext = dlg.GetWildcard().split('|')[2*dlg.GetFilterIndex()+1].strip('*')
78    if ext not in file:
79        file += ext
80    return file
81       
82def GetPowderPeaks(fileName):
83    'Read powder peaks from a file'
84    sind = lambda x: math.sin(x*math.pi/180.)
85    asind = lambda x: 180.*math.asin(x)/math.pi
86    Cuka = 1.54052
87    File = open(fileName,'Ur')
88    Comments = []
89    peaks = []
90    S = File.readline()
91    while S:
92        if S[:1] == '#':
93            Comments.append(S[:-1])
94        else:
95            item = S.split()
96            if len(item) == 1:
97                peaks.append([float(item[0]),1.0])
98            elif len(item) > 1:
99                peaks.append([float(item[0]),float(item[0])])
100        S = File.readline()
101    File.close()
102    if Comments:
103       print 'Comments on file:'
104       for Comment in Comments: print Comment
105    Peaks = []
106    if peaks[0][0] > peaks[-1][0]:          # d-spacings - assume CuKa
107        for peak in peaks:
108            dsp = peak[0]
109            sth = Cuka/(2.0*dsp)
110            if sth < 1.0:
111                tth = 2.0*asind(sth)
112            else:
113                break
114            Peaks.append([tth,peak[1],True,False,0,0,0,dsp,0.0])
115    else:                                   #2-thetas - assume Cuka (for now)
116        for peak in peaks:
117            tth = peak[0]
118            dsp = Cuka/(2.0*sind(tth/2.0))
119            Peaks.append([tth,peak[1],True,False,0,0,0,dsp,0.0])
120    return Comments,Peaks
121
122def CheckImageFile(G2frame,imagefile):
123    '''Get an new image file name if the specified one does not
124    exist
125
126    :param wx.Frame G2frame: main GSAS-II Frame and data object
127
128    :param str imagefile: name of image file
129
130    :returns: imagefile, if it exists, or the name of a file
131      that does exist or False if the user presses Cancel
132
133    '''
134    if not ospath.exists(imagefile):
135        dlg = wx.FileDialog(G2frame, 'Bad image file name; choose name', '.', '',\
136        'Any image file (*.edf;*.tif;*.tiff;*.mar*;*.avg;*.sum;*.img)\
137        |*.edf;*.tif;*.tiff;*.mar*;*.avg;*.sum;*.img|\
138        European detector file (*.edf)|*.edf|\
139        Any detector tif (*.tif;*.tiff)|*.tif;*.tiff|\
140        MAR file (*.mar*)|*.mar*|\
141        GE Image (*.avg;*.sum)|*.avg;*.sum|\
142        ADSC Image (*.img)|*.img|\
143        All files (*.*)|*.*',wx.OPEN|wx.CHANGE_DIR)
144        try:
145            dlg.SetFilename(''+ospath.split(imagefile)[1])
146            if dlg.ShowModal() == wx.ID_OK:
147                imagefile = dlg.GetPath()
148            else:
149                imagefile = False
150        finally:
151            dlg.Destroy()
152    return imagefile
153       
154def GetImageData(G2frame,imagefile,imageOnly=False):
155    '''Read an image with the file reader keyed by the
156    file extension
157
158    :param wx.Frame G2frame: main GSAS-II Frame and data object. Note: not used!
159
160    :param str imagefile: name of image file
161
162    :param bool imageOnly: If True return only the image,
163      otherwise  (default) return more (see below)
164
165    :returns: an image as a numpy array or a list of four items:
166      Comments, Data, Npix and the Image, as selected by imageOnly
167
168    '''
169    ext = ospath.splitext(imagefile)[1]
170    Comments = []
171    if ext == '.tif' or ext == '.tiff':
172        Comments,Data,Npix,Image = GetTifData(imagefile)
173    elif ext == '.edf':
174        Comments,Data,Npix,Image = GetEdfData(imagefile)
175    elif ext == '.img':
176        Comments,Data,Npix,Image = GetImgData(imagefile)
177        Image[0][0] = 0
178    elif ext == '.mar3450' or ext == '.mar2300':
179        Comments,Data,Npix,Image = GetMAR345Data(imagefile)
180    elif ext in ['.sum','.avg','']:
181        Comments,Data,Npix,Image = GetGEsumData(imagefile)
182    elif ext == '.G2img':
183        Comments,Data,Npix,Image = GetG2Image(imagefile)
184    if imageOnly:
185        return Image
186    else:
187        return Comments,Data,Npix,Image
188       
189def PutG2Image(filename,Comments,Data,Npix,image):
190    'Write an image as a python pickle'
191    File = open(filename,'wb')
192    cPickle.dump([Comments,Data,Npix,image],File,1)
193    File.close()
194    return
195   
196def GetG2Image(filename):
197    'Read an image as a python pickle'
198    File = open(filename,'rb')
199    Comments,Data,Npix,image = cPickle.load(File)
200    File.close()
201    return Comments,Data,Npix,image
202   
203def GetEdfData(filename,imageOnly=False):   
204    'Read European detector data edf file'
205    import struct as st
206    import array as ar
207    if not imageOnly:
208        print 'Read European detector data edf file: ',filename
209    File = open(filename,'rb')
210    fileSize = os.stat(filename).st_size
211    head = File.read(3072)
212    lines = head.split('\n')
213    sizexy = [0,0]
214    pixSize = [0,0]
215    cent = [0,0]
216    head = ['European detector data',]
217    for line in lines:
218        fields = line.split()
219        if 'Dim_1' in line:
220            sizexy[0] = int(fields[2])
221        elif 'Dim_2' in line:
222            sizexy[1] = int(fields[2])
223        elif 'DataType' in line:
224            dType = fields[2]
225        elif 'refined_wavelength' in line:
226            wave = float(fields[2])
227        elif 'Size' in line:
228            imSize = int(fields[2])
229        elif 'DataType' in lines:
230            dType = fields[2]
231        elif 'pixel_size_x' in line:
232            pixSize[0] = float(fields[2])
233        elif 'pixel_size_y' in line:
234            pixSize[1] = float(fields[2])
235        elif 'beam_center_x' in line:
236            cent[0] = float(fields[2])
237        elif 'beam_center_y' in line:
238            cent[1] = float(fields[2])
239        elif 'refined_distance' in line:
240            dist = float(fields[2])
241        if line:
242            head.append(line)
243    File.seek(fileSize-imSize)
244    if dType == 'UnsignedShort':       
245        image = np.array(ar.array('H',File.read(imSize)),dtype=np.int32)
246    image = np.reshape(image,(sizexy[1],sizexy[0]))
247    data = {'pixelSize':pixSize,'wavelength':wave,'distance':dist,'center':cent,'size':sizexy}
248    Npix = sizexy[0]*sizexy[1]
249    File.close()   
250    if imageOnly:
251        return image
252    else:
253        return head,data,Npix,image
254       
255def GetGEsumData(filename,imageOnly=False):
256    'Read SUM file as produced at 1-ID from G.E. images'
257    import struct as st
258    import array as ar
259    if not imageOnly:
260        print 'Read GE sum file: ',filename   
261    File = open(filename,'rb')
262    if '.sum' in filename:
263        head = ['GE detector sum data from APS 1-ID',]
264        sizexy = [2048,2048]
265    elif '.avg' in filename:
266        head = ['GE detector avg data from APS 1-ID',]
267        sizexy = [2048,2048]
268    else:
269        head = ['GE detector raw data from APS 1-ID',]
270        File.seek(18)
271        size,nframes = st.unpack('<ih',File.read(6))
272        sizexy = [2048,2048]
273        pos = 8192
274        File.seek(pos)
275    Npix = sizexy[0]*sizexy[1]
276    if '.sum' in filename:
277        image = np.array(ar.array('f',File.read(4*Npix)),dtype=np.int32)
278    elif '.avg' in filename:
279        image = np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)
280    else:
281        image = np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)
282        while nframes > 1:
283            image += np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)
284            nframes -= 1
285    image = np.reshape(image,(sizexy[1],sizexy[0]))
286    data = {'pixelSize':(200,200),'wavelength':0.15,'distance':250.0,'center':[204.8,204.8],'size':sizexy} 
287    File.close()   
288    if imageOnly:
289        return image
290    else:
291        return head,data,Npix,image
292       
293def GetImgData(filename,imageOnly=False):
294    'Read an ADSC image file'
295    import struct as st
296    import array as ar
297    if not imageOnly:
298        print 'Read ADSC img file: ',filename
299    File = open(filename,'rb')
300    head = File.read(511)
301    lines = head.split('\n')
302    head = []
303    center = [0,0]
304    for line in lines[1:-2]:
305        line = line.strip()[:-1]
306        if line:
307            if 'SIZE1' in line:
308                size = int(line.split('=')[1])
309                Npix = size*size
310            elif 'WAVELENGTH' in line:
311                wave = float(line.split('=')[1])
312            elif 'BIN' in line:
313                if line.split('=')[1] == '2x2':
314                    pixel=(102,102)
315                else:
316                    pixel = (51,51)
317            elif 'DISTANCE' in line:
318                distance = float(line.split('=')[1])
319            elif 'CENTER_X' in line:
320                center[0] = float(line.split('=')[1])
321            elif 'CENTER_Y' in line:
322                center[1] = float(line.split('=')[1])
323            head.append(line)
324    data = {'pixelSize':pixel,'wavelength':wave,'distance':distance,'center':center,'size':[size,size]}
325    image = []
326    row = 0
327    pos = 512
328    File.seek(pos)
329    image = np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)
330    image = np.reshape(image,(sizexy[1],sizexy[0]))
331#    image = np.zeros(shape=(size,size),dtype=np.int32)   
332#    while row < size:
333#        File.seek(pos)
334#        line = ar.array('H',File.read(2*size))
335#        image[row] = np.asarray(line)
336#        row += 1
337#        pos += 2*size
338    File.close()
339    if imageOnly:
340        return image
341    else:
342        return lines[1:-2],data,Npix,image
343       
344def GetMAR345Data(filename,imageOnly=False):
345    'Read a MAR-345 image plate image'
346    import array as ar
347    import struct as st
348    try:
349        import pack_f as pf
350    except:
351        msg = wx.MessageDialog(None, message="Unable to load the GSAS MAR image decompression, pack_f",
352                               caption="Import Error",
353                               style=wx.ICON_ERROR | wx.OK | wx.STAY_ON_TOP)
354        msg.ShowModal()
355        return None,None,None,None
356
357    if not imageOnly:
358        print 'Read Mar345 file: ',filename
359    File = open(filename,'rb')
360    head = File.read(4095)
361    numbers = st.unpack('<iiiiiiiiii',head[:40])
362    lines = head[128:].split('\n')
363    head = []
364    for line in lines:
365        line = line.strip()
366        if 'PIXEL' in line:
367            values = line.split()
368            pixel = (int(values[2]),int(values[4]))     #in microns
369        elif 'WAVELENGTH' in line:
370            wave = float(line.split()[1])
371        elif 'DISTANCE' in line:
372            distance = float(line.split()[1])           #in mm
373        elif 'CENTER' in line:
374            values = line.split()
375            center = [float(values[2])/10.,float(values[4])/10.]    #make in mm from pixels
376        if line: 
377            head.append(line)
378    data = {'pixelSize':pixel,'wavelength':wave,'distance':distance,'center':center}
379    for line in head:
380        if 'FORMAT' in line[0:6]:
381            items = line.split()
382            size = int(items[1])
383            Npix = size*size
384    pos = 4096
385    data['size'] = [size,size]
386    File.seek(pos)
387    line = File.read(8)
388    while 'CCP4' not in line:       #get past overflow list for now
389        line = File.read(8)
390        pos += 8
391    pos += 37
392    File.seek(pos)
393    raw = File.read()
394    File.close()
395    image = np.zeros(shape=(size,size),dtype=np.int32)
396    image = np.flipud(pf.pack_f(len(raw),raw,size,image).T)  #transpose to get it right way around & flip
397    if imageOnly:
398        return image
399    else:
400        return head,data,Npix,image
401
402def GetTifData(filename,imageOnly=False):
403    '''Read an image in a pseudo-tif format,
404    as produced by a wide variety of software, almost always
405    incorrectly in some way.
406    '''
407    import struct as st
408    import array as ar
409    import ReadMarCCDFrame as rmf
410    File = open(filename,'rb')
411    dataType = 5
412    try:
413        Meta = open(filename+'.metadata','Ur')
414        head = Meta.readlines()
415        for line in head:
416            line = line.strip()
417            if 'dataType=' in line:
418                dataType = int(line.split('=')[1])
419        Meta.close()
420    except IOError:
421        print 'no metadata file found - will try to read file anyway'
422        head = ['no metadata file found',]
423       
424    tag = File.read(2)
425    byteOrd = '<'
426    if tag == 'II' and int(st.unpack('<h',File.read(2))[0]) == 42:     #little endian
427        IFD = int(st.unpack(byteOrd+'i',File.read(4))[0])
428    elif tag == 'MM' and int(st.unpack('>h',File.read(2))[0]) == 42:   #big endian
429        byteOrd = '>'
430        IFD = int(st.unpack(byteOrd+'i',File.read(4))[0])       
431    else:
432        lines = ['not a detector tiff file',]
433        return lines,0,0,0
434    File.seek(IFD)                                                  #get number of directory entries
435    NED = int(st.unpack(byteOrd+'h',File.read(2))[0])
436    IFD = {}
437    for ied in range(NED):
438        Tag,Type = st.unpack(byteOrd+'Hh',File.read(4))
439        nVal = st.unpack(byteOrd+'i',File.read(4))[0]
440        if Type == 1:
441            Value = st.unpack(byteOrd+nVal*'b',File.read(nVal))
442        elif Type == 2:
443            Value = st.unpack(byteOrd+'i',File.read(4))
444        elif Type == 3:
445            Value = st.unpack(byteOrd+nVal*'h',File.read(nVal*2))
446            x = st.unpack(byteOrd+nVal*'h',File.read(nVal*2))
447        elif Type == 4:
448            Value = st.unpack(byteOrd+nVal*'i',File.read(nVal*4))
449        elif Type == 5:
450            Value = st.unpack(byteOrd+nVal*'i',File.read(nVal*4))
451        elif Type == 11:
452            Value = st.unpack(byteOrd+nVal*'f',File.read(nVal*4))
453        IFD[Tag] = [Type,nVal,Value]
454        #print Tag,IFD[Tag]
455    sizexy = [IFD[256][2][0],IFD[257][2][0]]
456    [nx,ny] = sizexy
457    Npix = nx*ny
458    if 34710 in IFD:
459        if not imageOnly:
460            print 'Read MAR CCD tiff file: ',filename
461        marFrame = rmf.marFrame(File,byteOrd,IFD)
462        image = np.flipud(np.array(np.asarray(marFrame.image),dtype=np.int32))
463        tifType = marFrame.filetitle
464        pixy = (marFrame.pixelsizeX/1000.0,marFrame.pixelsizeY/1000.0)
465        head = marFrame.outputHead()
466    elif 272 in IFD:
467        ifd = IFD[272]
468        File.seek(ifd[2][0])
469        S = File.read(ifd[1])
470        if 'PILATUS' in S:
471            tifType = 'Pilatus'
472            dataType = 0
473            pixy = (172,172)
474            File.seek(4096)
475            if not imageOnly:
476                print 'Read Pilatus tiff file: ',filename
477            image = ar.array('L',File.read(4*Npix))
478            image = np.array(np.asarray(image),dtype=np.int32)
479        else:
480            if IFD[258][2][0] == 16:
481                print sizexy
482                tifType = 'GE'
483                pixy = (200,200)
484                File.seek(8)
485                if not imageOnly:
486                    print 'Read GE-detector tiff file: ',filename
487                image = np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)
488           
489    elif 262 in IFD and IFD[262][2][0] > 4:
490        tifType = 'DND'
491        pixy = (158,158)
492        File.seek(512)
493        if not imageOnly:
494            print 'Read DND SAX/WAX-detector tiff file: ',filename
495        image = np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)
496    elif sizexy == [1536,1536]:
497        tifType = 'APS Gold'
498        pixy = (150,150)
499        File.seek(64)
500        if not imageOnly:
501            print 'Read Gold tiff file:',filename
502        image = np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)
503    elif sizexy == [2048,2048] or sizexy == [1024,1024] or sizexy == [3072,3072]:
504        if IFD[273][2][0] == 8:
505            if IFD[258][2][0] == 32:
506                tifType = 'PE'
507                pixy = (200,200)
508                File.seek(8)
509                if not imageOnly:
510                    print 'Read APS PE-detector tiff file: ',filename
511                if dataType == 5:
512                    image = np.array(ar.array('f',File.read(4*Npix)),dtype=np.float32)
513                else:
514                    image = np.array(ar.array('I',File.read(4*Npix)),dtype=np.int32)
515               
516        elif IFD[273][2][0] == 4096:
517            if sizexy[0] == 3072:
518                pixy =  (73,73)
519                tifType = 'MAR225'           
520            else:
521                pixy = (158,158)
522                tifType = 'MAR325'           
523            File.seek(4096)
524            if not imageOnly:
525                print 'Read MAR CCD tiff file: ',filename
526            image = np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)
527        elif IFD[273][2][0] == 512:
528            tiftype = '11-ID-C'
529            pixy = [200,200]
530            File.seek(512)
531            if not imageOnly:
532                print 'Read 11-ID-C tiff file: ',filename
533            image = np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)           
534    elif sizexy == [4096,4096]:
535        if IFD[273][2][0] == 8:
536            if IFD[258][2][0] == 16:
537                tifType = 'scanCCD'
538                pixy = (9,9)
539                File.seek(8)
540                if not imageOnly:
541                    print 'Read APS scanCCD tiff file: ',filename
542                image = np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)
543        elif IFD[273][2][0] == 4096:
544            tifType = 'Rayonix'
545            pixy = (73.242,73.242)
546            File.seek(4096)
547            if not imageOnly:
548                print 'Read Rayonix MX300HE tiff file: ',filename
549            image = np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)
550#    elif sizexy == [960,960]:
551#        tiftype = 'PE-BE'
552#        pixy = (200,200)
553#        File.seek(8)
554#        if not imageOnly:
555#            print 'Read Gold tiff file:',filename
556#        image = np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)
557           
558    else:
559        lines = ['not a known detector tiff file',]
560        return lines,0,0,0
561       
562    image = np.reshape(image,(sizexy[1],sizexy[0]))
563    center = [pixy[0]*sizexy[0]/2000,pixy[1]*sizexy[1]/2000]
564    data = {'pixelSize':pixy,'wavelength':0.10,'distance':100.0,'center':center,'size':sizexy}
565    File.close()   
566    if imageOnly:
567        return image
568    else:
569        return head,data,Npix,image
570   
571#def GetTifData(filename,imageOnly=False):
572#    import struct as st
573#    import array as ar
574#    File = open(filename,'rb')
575#    dataType = 5
576#    try:
577#        Meta = open(filename+'.metadata','Ur')
578#        head = Meta.readlines()
579#        for line in head:
580#            line = line.strip()
581#            if 'dataType=' in line:
582#                dataType = int(line.split('=')[1])
583#        Meta.close()
584#    except IOError:
585#        print 'no metadata file found - will try to read file anyway'
586#        head = ['no metadata file found',]
587#       
588#    tag = File.read(2)
589#    byteOrd = '<'
590#    if tag == 'II' and int(st.unpack('<h',File.read(2))[0]) == 42:     #little endian
591#        IFD = int(st.unpack(byteOrd+'i',File.read(4))[0])
592#    elif tag == 'MM' and int(st.unpack('>h',File.read(2))[0]) == 42:   #big endian
593#        byteOrd = '>'
594#        IFD = int(st.unpack(byteOrd+'i',File.read(4))[0])       
595#    else:
596#        lines = ['not a detector tiff file',]
597#        return lines,0,0,0
598#    File.seek(IFD)                                                  #get number of directory entries
599#    NED = int(st.unpack(byteOrd+'h',File.read(2))[0])
600#    IFD = {}
601#    for ied in range(NED):
602#        Tag,Type = st.unpack(byteOrd+'Hh',File.read(4))
603#        nVal = st.unpack(byteOrd+'i',File.read(4))[0]
604#        if Type == 1:
605#            Value = st.unpack(byteOrd+nVal*'b',File.read(nVal))
606#        elif Type == 2:
607#            Value = st.unpack(byteOrd+'i',File.read(4))
608#        elif Type == 3:
609#            Value = st.unpack(byteOrd+nVal*'h',File.read(nVal*2))
610#            x = st.unpack(byteOrd+nVal*'h',File.read(nVal*2))
611#        elif Type == 4:
612#            Value = st.unpack(byteOrd+nVal*'i',File.read(nVal*4))
613#        elif Type == 5:
614#            Value = st.unpack(byteOrd+nVal*'i',File.read(nVal*4))
615#        elif Type == 11:
616#            Value = st.unpack(byteOrd+nVal*'f',File.read(nVal*4))
617#        IFD[Tag] = [Type,nVal,Value]
618##        print Tag,IFD[Tag]
619#    sizexy = [IFD[256][2][0],IFD[257][2][0]]
620#    [nx,ny] = sizexy
621#    Npix = nx*ny
622#    if 272 in IFD:
623#        ifd = IFD[272]
624#        File.seek(ifd[2][0])
625#        S = File.read(ifd[1])
626#        if 'PILATUS' in S:
627#            tifType = 'Pilatus'
628#            dataType = 0
629#            pixy = (172,172)
630#            File.seek(4096)
631#            if not imageOnly:
632#                print 'Read Pilatus tiff file: ',filename
633#            image = ar.array('L',File.read(4*Npix))
634#            image = np.array(np.asarray(image),dtype=np.int32)
635#    elif 262 in IFD and IFD[262][2][0] > 4:
636#        tifType = 'DND'
637#        pixy = (158,158)
638#        File.seek(512)
639#        if not imageOnly:
640#            print 'Read DND SAX/WAX-detector tiff file: ',filename
641#        image = np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)
642#    elif sizexy == [1536,1536]:
643#        tifType = 'APS Gold'
644#        pixy = (150,150)
645#        File.seek(64)
646#        if not imageOnly:
647#            print 'Read Gold tiff file:',filename
648#        image = np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)
649#    elif sizexy == [2048,2048] or sizexy == [1024,1024] or sizexy == [3072,3072]:
650#        if IFD[273][2][0] == 8:
651#            if IFD[258][2][0] == 32:
652#                tifType = 'PE'
653#                pixy = (200,200)
654#                File.seek(8)
655#                if not imageOnly:
656#                    print 'Read APS PE-detector tiff file: ',filename
657#                if dataType == 5:
658#                    image = np.array(ar.array('f',File.read(4*Npix)),dtype=np.float32)
659#                else:
660#                    image = np.array(ar.array('I',File.read(4*Npix)),dtype=np.int32)
661#        elif IFD[273][2][0] == 4096:
662#            if sizexy[0] == 3072:
663#                pixy =  (73,73)
664#                tifType = 'MAR225'           
665#            else:
666#                pixy = (158,158)
667#                tifType = 'MAR325'           
668#            File.seek(4096)
669#            if not imageOnly:
670#                print 'Read MAR CCD tiff file: ',filename
671#            image = np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)
672#        elif IFD[273][2][0] == 512:
673#            tiftype = '11-ID-C'
674#            pixy = [200,200]
675#            File.seek(512)
676#            if not imageOnly:
677#                print 'Read 11-ID-C tiff file: ',filename
678#            image = np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)           
679#    elif sizexy == [4096,4096]:
680#        if IFD[273][2][0] == 8:
681#            if IFD[258][2][0] == 16:
682#                tifType = 'scanCCD'
683#                pixy = (9,9)
684#                File.seek(8)
685#                if not imageOnly:
686#                    print 'Read APS scanCCD tiff file: ',filename
687#                image = np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)
688#        elif IFD[273][2][0] == 4096:
689#            tifType = 'Rayonix'
690#            pixy = (73.242,73.242)
691#            File.seek(4096)
692#            if not imageOnly:
693#                print 'Read Rayonix MX300HE tiff file: ',filename
694#            image = np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)
695##    elif sizexy == [960,960]:
696##        tiftype = 'PE-BE'
697##        pixy = (200,200)
698##        File.seek(8)
699##        if not imageOnly:
700##            print 'Read Gold tiff file:',filename
701##        image = np.array(ar.array('H',File.read(2*Npix)),dtype=np.int32)
702#           
703#    else:
704#        lines = ['not a known detector tiff file',]
705#        return lines,0,0,0
706#       
707#    image = np.reshape(image,(sizexy[1],sizexy[0]))
708#    center = [pixy[0]*sizexy[0]/2000,pixy[1]*sizexy[1]/2000]
709#    data = {'pixelSize':pixy,'wavelength':0.10,'distance':100.0,'center':center,'size':sizexy}
710#    File.close()   
711#    if imageOnly:
712#        return image
713#    else:
714#        return head,data,Npix,image
715#   
716def ProjFileOpen(G2frame):
717    'Read a GSAS-II project file and load into the G2 data tree'
718    file = open(G2frame.GSASprojectfile,'rb')
719    print 'load from file: ',G2frame.GSASprojectfile
720    G2frame.SetTitle("GSAS-II data tree: "+
721                     os.path.split(G2frame.GSASprojectfile)[1])
722    wx.BeginBusyCursor()
723    try:
724        while True:
725            try:
726                data = cPickle.load(file)
727            except EOFError:
728                break
729            datum = data[0]
730           
731            Id = G2frame.PatternTree.AppendItem(parent=G2frame.root,text=datum[0])
732            if 'PWDR' in datum[0]:               
733                if 'ranId' not in datum[1][0]: # patch: add random Id if not present
734                    datum[1][0]['ranId'] = ran.randint(0,sys.maxint)
735                G2frame.PatternTree.SetItemPyData(Id,datum[1][:3])  #temp. trim off junk (patch?)
736            elif datum[0].startswith('HKLF'): 
737                if 'ranId' not in datum[1][0]: # patch: add random Id if not present
738                    datum[1][0]['ranId'] = ran.randint(0,sys.maxint)
739                G2frame.PatternTree.SetItemPyData(Id,datum[1])
740            else:
741                G2frame.PatternTree.SetItemPyData(Id,datum[1])
742            for datus in data[1:]:
743                sub = G2frame.PatternTree.AppendItem(Id,datus[0])
744#patch
745                if datus[0] == 'Instrument Parameters' and len(datus[1]) == 1:
746                    if 'PWDR' in datum[0]:
747                        datus[1] = [dict(zip(datus[1][3],zip(datus[1][0],datus[1][1],datus[1][2]))),{}]
748                    else:
749                        datus[1] = [dict(zip(datus[1][2],zip(datus[1][0],datus[1][1]))),{}]
750                    for item in datus[1][0]:               #zip makes tuples - now make lists!
751                        datus[1][0][item] = list(datus[1][0][item])
752#end patch
753                G2frame.PatternTree.SetItemPyData(sub,datus[1])
754            if 'IMG' in datum[0]:                   #retrieve image default flag & data if set
755                Data = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,Id,'Image Controls'))
756                if Data['setDefault']:
757                    G2frame.imageDefault = Data               
758        file.close()
759        print('project load successful')
760        G2frame.NewPlot = True
761    except:
762        msg = wx.MessageDialog(G2frame,message="Error reading file "+
763            str(G2frame.GSASprojectfile)+". This is not a GSAS-II .gpx file",
764            caption="Load Error",style=wx.ICON_ERROR | wx.OK | wx.STAY_ON_TOP)
765        msg.ShowModal()
766    finally:
767        wx.EndBusyCursor()
768   
769def ProjFileSave(G2frame):
770    'Save a GSAS-II project file'
771    if not G2frame.PatternTree.IsEmpty():
772        file = open(G2frame.GSASprojectfile,'wb')
773        print 'save to file: ',G2frame.GSASprojectfile
774        wx.BeginBusyCursor()
775        try:
776            item, cookie = G2frame.PatternTree.GetFirstChild(G2frame.root)
777            while item:
778                data = []
779                name = G2frame.PatternTree.GetItemText(item)
780                data.append([name,G2frame.PatternTree.GetItemPyData(item)])
781                item2, cookie2 = G2frame.PatternTree.GetFirstChild(item)
782                while item2:
783                    name = G2frame.PatternTree.GetItemText(item2)
784                    data.append([name,G2frame.PatternTree.GetItemPyData(item2)])
785                    item2, cookie2 = G2frame.PatternTree.GetNextChild(item, cookie2)                           
786                item, cookie = G2frame.PatternTree.GetNextChild(G2frame.root, cookie)                           
787                cPickle.dump(data,file,1)
788            file.close()
789        finally:
790            wx.EndBusyCursor()
791        print('project save successful')
792
793def SaveIntegration(G2frame,PickId,data):
794    'Save image integration results as powder pattern(s)'
795    azms = G2frame.Integrate[1]
796    X = G2frame.Integrate[2][:-1]
797    Xminmax = [X[0],X[-1]]
798    N = len(X)
799    Id = G2frame.PatternTree.GetItemParent(PickId)
800    name = G2frame.PatternTree.GetItemText(Id)
801    name = name.replace('IMG ','PWDR ')
802    Comments = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,Id, 'Comments'))
803    names = ['Type','Lam','Zero','Polariz.','U','V','W','X','Y','SH/L','Azimuth'] 
804    codes = [0 for i in range(12)]
805    LRazm = data['LRazimuth']
806    Azms = []
807    if data['fullIntegrate'] and data['outAzimuths'] == 1:
808        Azms = [45.0,]                              #a poor man's average?
809    else:
810        for i,azm in enumerate(azms[:-1]):
811            Azms.append((azms[i+1]+azm)/2.)
812    for i,azm in enumerate(azms[:-1]):
813        item, cookie = G2frame.PatternTree.GetFirstChild(G2frame.root)
814        Id = 0
815        while item:
816            Name = G2frame.PatternTree.GetItemText(item)
817            if name == Name:
818                Id = item
819            item, cookie = G2frame.PatternTree.GetNextChild(G2frame.root, cookie)
820        parms = ['PXC',data['wavelength'],0.0,0.99,1.0,-0.10,0.4,0.30,1.0,0.0001,Azms[i]]    #set polarization for synchrotron radiation!
821        Y = G2frame.Integrate[0][i]
822        W = 1./Y                    #probably not true
823        Sample = G2pdG.SetDefaultSample()
824        Sample['Gonio. radius'] = data['distance']
825        Sample['Omega'] = data['GonioAngles'][0]
826        Sample['Chi'] = data['GonioAngles'][1]
827        Sample['Phi'] = data['GonioAngles'][2]
828        if Id:
829            G2frame.PatternTree.SetItemPyData(G2gd.GetPatternTreeItemId(G2frame,Id, 'Comments'),Comments)                   
830            G2frame.PatternTree.SetItemPyData(G2gd.GetPatternTreeItemId(G2frame,Id,'Limits'),[tuple(Xminmax),Xminmax])
831            G2frame.PatternTree.SetItemPyData(G2gd.GetPatternTreeItemId(G2frame,Id,'Background'),[['chebyschev',1,3,1.0,0.0,0.0],
832                            {'nDebye':0,'debyeTerms':[],'nPeaks':0,'peaksList':[]}])
833            inst = [dict(zip(names,zip(parms,parms,codes))),{}]
834            for item in inst[0]:
835                inst[0][item] = list(inst[0][item])
836            G2frame.PatternTree.SetItemPyData(G2gd.GetPatternTreeItemId(G2frame,Id,'Instrument Parameters'),inst)
837            G2frame.PatternTree.SetItemPyData(G2gd.GetPatternTreeItemId(G2frame,Id,'Peak List'),[])
838            G2frame.PatternTree.SetItemPyData(G2gd.GetPatternTreeItemId(G2frame,Id,'Index Peak List'),[])
839            G2frame.PatternTree.SetItemPyData(G2gd.GetPatternTreeItemId(G2frame,Id,'Unit Cells List'),[])             
840            G2frame.PatternTree.SetItemPyData(G2gd.GetPatternTreeItemId(G2frame,Id,'Reflection Lists'),{})             
841        else:
842            Id = G2frame.PatternTree.AppendItem(parent=G2frame.root,text=name+" Azm= %.2f"%(Azms[i]))
843            G2frame.PatternTree.SetItemPyData(G2frame.PatternTree.AppendItem(Id,text='Comments'),Comments)                   
844            G2frame.PatternTree.SetItemPyData(G2frame.PatternTree.AppendItem(Id,text='Limits'),[tuple(Xminmax),Xminmax])
845            G2frame.PatternTree.SetItemPyData(G2frame.PatternTree.AppendItem(Id,text='Background'),[['chebyschev',1,3,1.0,0.0,0.0],
846                            {'nDebye':0,'debyeTerms':[],'nPeaks':0,'peaksList':[]}])
847            inst = [dict(zip(names,zip(parms,parms,codes))),{}]
848            for item in inst[0]:
849                inst[0][item] = list(inst[0][item])
850            G2frame.PatternTree.SetItemPyData(G2frame.PatternTree.AppendItem(Id,text='Instrument Parameters'),inst)
851            G2frame.PatternTree.SetItemPyData(G2frame.PatternTree.AppendItem(Id,text='Sample Parameters'),Sample)
852            G2frame.PatternTree.SetItemPyData(G2frame.PatternTree.AppendItem(Id,text='Peak List'),[])
853            G2frame.PatternTree.SetItemPyData(G2frame.PatternTree.AppendItem(Id,text='Index Peak List'),[])
854            G2frame.PatternTree.SetItemPyData(G2frame.PatternTree.AppendItem(Id,text='Unit Cells List'),[])
855            G2frame.PatternTree.SetItemPyData(G2frame.PatternTree.AppendItem(Id,text='Reflection Lists'),{})
856            valuesdict = {
857                'wtFactor':1.0,
858                'Dummy':False,
859                'ranId':ran.randint(0,sys.maxint),
860                }
861        G2frame.PatternTree.SetItemPyData(
862            Id,[valuesdict,
863                [np.array(X),np.array(Y),np.array(W),np.zeros(N),np.zeros(N),np.zeros(N)]])
864    G2frame.PatternTree.SelectItem(Id)
865    G2frame.PatternTree.Expand(Id)
866    G2frame.PatternId = Id
867           
868def powderFxyeSave(G2frame,exports,powderfile):
869    'Save a powder histogram as a GSAS FXYE file'
870    head,tail = ospath.split(powderfile)
871    name,ext = tail.split('.')
872    for i,export in enumerate(exports):
873        filename = ospath.join(head,name+'-%03d.'%(i)+ext)
874        prmname = filename.strip(ext)+'prm'
875        prm = open(prmname,'w')      #old style GSAS parm file
876        PickId = G2gd.GetPatternTreeItemId(G2frame, G2frame.root, export)
877        Inst = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame, \
878            PickId, 'Instrument Parameters'))[0]
879        prm.write( '            123456789012345678901234567890123456789012345678901234567890        '+'\n')
880        prm.write( 'INS   BANK      1                                                               '+'\n')
881        prm.write(('INS   HTYPE   %sR                                                              '+'\n')%(Inst['Type'][0]))
882        if 'Lam1' in Inst:              #Ka1 & Ka2
883            prm.write(('INS  1 ICONS%10.7f%10.7f    0.0000               0.990    0     0.500   '+'\n')%(Inst['Lam1'][0],Inst['Lam2'][0]))
884        elif 'Lam' in Inst:             #single wavelength
885            prm.write(('INS  1 ICONS%10.7f%10.7f    0.0000               0.990    0     0.500   '+'\n')%(Inst['Lam'][1],0.0))
886        prm.write( 'INS  1 IRAD     0                                                               '+'\n')
887        prm.write( 'INS  1I HEAD                                                                    '+'\n')
888        prm.write( 'INS  1I ITYP    0    0.0000  180.0000         1                                 '+'\n')
889        prm.write(('INS  1DETAZM%10.3f                                                          '+'\n')%(Inst['Azimuth'][0]))
890        prm.write( 'INS  1PRCF1     3    8   0.00100                                                '+'\n')
891        prm.write(('INS  1PRCF11     %15.6g%15.6g%15.6g%15.6g   '+'\n')%(Inst['U'][1],Inst['V'][1],Inst['W'][1],0.0))
892        prm.write(('INS  1PRCF12     %15.6g%15.6g%15.6g%15.6g   '+'\n')%(Inst['X'][1],Inst['Y'][1],Inst['SH/L'][1]/2.,Inst['SH/L'][1]/2.))
893        prm.close()
894        file = open(filename,'w')
895        print 'save powder pattern to file: ',filename
896        x,y,w,yc,yb,yd = G2frame.PatternTree.GetItemPyData(PickId)[1]
897        file.write(powderfile+'\n')
898        file.write('Instrument parameter file:'+ospath.split(prmname)[1]+'\n')
899        file.write('BANK 1 %d %d CONS %.2f %.2f 0 0 FXYE\n'%(len(x),len(x),\
900            100.*x[0],100.*(x[1]-x[0])))
901        s = list(np.sqrt(1./np.array(w)))       
902        XYW = zip(x,y,s)
903        for X,Y,S in XYW:
904            file.write("%15.6g %15.6g %15.6g\n" % (100.*X,Y,max(S,1.0)))
905        file.close()
906        print 'powder pattern file '+filename+' written'
907       
908def powderXyeSave(G2frame,exports,powderfile):
909    'Save a powder histogram as a Topas XYE file'
910    head,tail = ospath.split(powderfile)
911    name,ext = tail.split('.')
912    for i,export in enumerate(exports):
913        filename = ospath.join(head,name+'-%03d.'%(i)+ext)
914        PickId = G2gd.GetPatternTreeItemId(G2frame, G2frame.root, export)
915        file = open(filename,'w')
916        file.write('#%s\n'%(export))
917        print 'save powder pattern to file: ',filename
918        x,y,w,yc,yb,yd = G2frame.PatternTree.GetItemPyData(PickId)[1]
919        s = list(np.sqrt(1./np.array(w)))       
920        XYW = zip(x,y,s)
921        for X,Y,W in XYW:
922            file.write("%15.6g %15.6g %15.6g\n" % (X,Y,W))
923        file.close()
924        print 'powder pattern file '+filename+' written'
925       
926def PDFSave(G2frame,exports):
927    'Save a PDF G(r) and S(Q) in column formats'
928    for export in exports:
929        PickId = G2gd.GetPatternTreeItemId(G2frame, G2frame.root, export)
930        SQname = 'S(Q)'+export[4:]
931        GRname = 'G(R)'+export[4:]
932        sqfilename = ospath.join(G2frame.dirname,export.replace(' ','_')[5:]+'.sq')
933        grfilename = ospath.join(G2frame.dirname,export.replace(' ','_')[5:]+'.gr')
934        sqId = G2gd.GetPatternTreeItemId(G2frame, PickId, SQname)
935        grId = G2gd.GetPatternTreeItemId(G2frame, PickId, GRname)
936        sqdata = np.array(G2frame.PatternTree.GetItemPyData(sqId)[1][:2]).T
937        grdata = np.array(G2frame.PatternTree.GetItemPyData(grId)[1][:2]).T
938        sqfile = open(sqfilename,'w')
939        grfile = open(grfilename,'w')
940        sqfile.write('#T S(Q) %s\n'%(export))
941        grfile.write('#T G(R) %s\n'%(export))
942        sqfile.write('#L Q     S(Q)\n')
943        grfile.write('#L R     G(R)\n')
944        for q,sq in sqdata:
945            sqfile.write("%15.6g %15.6g\n" % (q,sq))
946        sqfile.close()
947        for r,gr in grdata:
948            grfile.write("%15.6g %15.6g\n" % (r,gr))
949        grfile.close()
950   
951def PeakListSave(G2frame,file,peaks):
952    'Save powder peaks to a data file'
953    print 'save peak list to file: ',G2frame.peaklistfile
954    if not peaks:
955        dlg = wx.MessageDialog(G2frame, 'No peaks!', 'Nothing to save!', wx.OK)
956        try:
957            result = dlg.ShowModal()
958        finally:
959            dlg.Destroy()
960        return
961    for peak in peaks:
962        file.write("%10.4f %12.2f %10.3f %10.3f \n" % \
963            (peak[0],peak[2],peak[4],peak[6]))
964    print 'peak list saved'
965             
966def IndexPeakListSave(G2frame,peaks):
967    'Save powder peaks from the indexing list'
968    file = open(G2frame.peaklistfile,'wa')
969    print 'save index peak list to file: ',G2frame.peaklistfile
970    wx.BeginBusyCursor()
971    try:
972        if not peaks:
973            dlg = wx.MessageDialog(G2frame, 'No peaks!', 'Nothing to save!', wx.OK)
974            try:
975                result = dlg.ShowModal()
976            finally:
977                dlg.Destroy()
978            return
979        for peak in peaks:
980            file.write("%12.6f\n" % (peak[7]))
981        file.close()
982    finally:
983        wx.EndBusyCursor()
984    print 'index peak list saved'
985   
986def SetNewPhase(Name='New Phase',SGData=None,cell=None):
987    '''Create a new phase with default values for various parameters
988
989    :param str Name: Name for new Phase
990
991    :param dict SGData: space group data from :func:`GSASIIspc:SpcGroup`;
992      defaults to data for P 1
993
994    :param list cell: unit cell parameter list; defaults to
995      [1.0,1.0,1.0,90.,90,90.,1.]
996
997    '''
998    if SGData is None: SGData = G2spc.SpcGroup('P 1')[1]
999    if cell is None: cell=[1.0,1.0,1.0,90.,90,90.,1.]
1000    phaseData = {
1001        'ranId':ran.randint(0,sys.maxint),
1002        'General':{
1003            'Name':Name,
1004            'Type':'nuclear',
1005            'SGData':SGData,
1006            'Cell':[False,]+cell,
1007            'Pawley dmin':1.0,
1008            'Data plot type':'None',
1009            'SH Texture':{
1010                'Order':0,
1011                'Model':'cylindrical',
1012                'Sample omega':[False,0.0],
1013                'Sample chi':[False,0.0],
1014                'Sample phi':[False,0.0],
1015                'SH Coeff':[False,{}],
1016                'SHShow':False,
1017                'PFhkl':[0,0,1],
1018                'PFxyz':[0,0,1],
1019                'PlotType':'Pole figure'}},
1020        'Atoms':[],
1021        'Drawing':{},
1022        'Histograms':{},
1023        'Pawley ref':[],
1024        'RBModels':{},
1025        }
1026    return phaseData
1027   
1028def ReadEXPPhase(G2frame,filename):
1029    '''Read a phase from a GSAS .EXP file.
1030    Called in the GSAS phase import routine (see imports/G2phase.py)
1031    '''
1032    shModels = ['cylindrical','none','shear - 2/m','rolling - mmm']
1033    textureData = {'Order':0,'Model':'cylindrical','Sample omega':[False,0.0],
1034        'Sample chi':[False,0.0],'Sample phi':[False,0.0],'SH Coeff':[False,{}],
1035        'SHShow':False,'PFhkl':[0,0,1],'PFxyz':[0,0,1],'PlotType':'Pole figure'}
1036    shNcof = 0
1037    file = open(filename, 'Ur')
1038    S = 1
1039    Expr = [{},{},{},{},{},{},{},{},{}]
1040    while S:
1041        S = file.readline()
1042        if 'EXPR NPHAS' in S[:12]:
1043            Num = S[12:-1].count('0')
1044            NPhas = S[12:-1].split()
1045        if 'CRS' in S[:3]:
1046            N = int(S[3:4])-1
1047            Expr[N][S[:12]] = S[12:-1]
1048    file.close()
1049    PNames = []
1050    for n,N in enumerate(NPhas):
1051        if N != '0':
1052            result = n
1053            key = 'CRS'+str(n+1)+'    PNAM'
1054            PNames.append(Expr[n][key])
1055    if Num < 8:
1056        dlg = wx.SingleChoiceDialog(G2frame, 'Which phase to read?', 'Read phase data', PNames, wx.CHOICEDLG_STYLE)
1057        try:
1058            if dlg.ShowModal() == wx.ID_OK:
1059                result = dlg.GetSelection()
1060        finally:
1061            dlg.Destroy()       
1062    EXPphase = Expr[result]
1063    keyList = EXPphase.keys()
1064    keyList.sort()
1065    SGData = {}
1066    if NPhas[result] == '1':
1067        Ptype = 'nuclear'
1068    elif NPhas[result] in ['2','3']:
1069        Ptype = 'magnetic'
1070    elif NPhas[result] == '4':
1071        Ptype = 'macromolecular'
1072    elif NPhas[result] == '10':
1073        Ptype = 'Pawley'
1074    for key in keyList:
1075        if 'PNAM' in key:
1076           PhaseName = EXPphase[key].strip()
1077        elif 'ABC   ' in key:
1078            abc = [float(EXPphase[key][:10]),float(EXPphase[key][10:20]),float(EXPphase[key][20:30])]                       
1079        elif 'ANGLES' in key:
1080            angles = [float(EXPphase[key][:10]),float(EXPphase[key][10:20]),float(EXPphase[key][20:30])]                                               
1081        elif 'SG SYM' in key:
1082            SpGrp = EXPphase[key][:15].strip()
1083            E,SGData = G2spc.SpcGroup(SpGrp)
1084            if E:
1085                E,SGData = G2spc.SpcGroup('P 1') # unlikely to need this
1086        elif 'OD    ' in key:
1087            SHdata = EXPphase[key].split() # may not have all 9 values
1088            SHvals = 9*[0]
1089            for i in range(9):
1090                try:
1091                    float(SHdata[i])
1092                    SHvals[i] = SHdata[i]
1093                except:
1094                    pass
1095            textureData['Order'] = int(SHvals[0])
1096            textureData['Model'] = shModels[int(SHvals[2])]
1097            textureData['Sample omega'] = [False,float(SHvals[6])]
1098            textureData['Sample chi'] = [False,float(SHvals[7])]
1099            textureData['Sample phi'] = [False,float(SHvals[8])]
1100            shNcof = int(SHvals[1])
1101    Atoms = []
1102    if Ptype == 'nuclear':
1103        for key in keyList:
1104            if 'AT' in key:
1105                if key[11:] == 'A':
1106                    S = EXPphase[key]
1107                elif key[11:] == 'B':
1108                    S += EXPphase[key]
1109                    Atom = [S[50:58].strip(),S[:10].strip().capitalize(),'',
1110                        float(S[10:20]),float(S[20:30]),float(S[30:40]),
1111                        float(S[40:50]),'',int(S[60:62]),S[130:131]]
1112                    if Atom[9] == 'I':
1113                        Atom += [float(S[68:78]),0.,0.,0.,0.,0.,0.]
1114                    elif Atom[9] == 'A':
1115                        Atom += [0.0,float(S[68:78]),float(S[78:88]),
1116                            float(S[88:98]),float(S[98:108]),
1117                            float(S[108:118]),float(S[118:128])]
1118                    XYZ = Atom[3:6]
1119                    Atom[7],Atom[8] = G2spc.SytSym(XYZ,SGData)
1120                    Atom.append(ran.randint(0,sys.maxint))
1121                    Atoms.append(Atom)
1122    elif Ptype == 'macromolecular':
1123        for key in keyList:
1124            if 'AT' in key[6:8]:
1125                S = EXPphase[key]
1126                Atom = [S[56:60],S[50:54].strip().upper(),S[54:56],
1127                    S[46:51].strip(),S[:8].strip().capitalize(),'',
1128                    float(S[16:24]),float(S[24:32]),float(S[32:40]),
1129                    float(S[8:16]),'1',1,'I',float(S[40:46]),0,0,0,0,0,0]
1130                XYZ = Atom[6:9]
1131                Atom[10],Atom[11] = G2spc.SytSym(XYZ,SGData)
1132                Atom.append(ran.randint(0,sys.maxint))
1133                Atoms.append(Atom)
1134    Volume = G2lat.calc_V(G2lat.cell2A(abc+angles))
1135    if shNcof:
1136        shCoef = {}
1137        nRec = [i+1 for i in range((shNcof-1)/6+1)]
1138        for irec in nRec:
1139            ODkey = keyList[0][:6]+'OD'+'%3dA'%(irec)
1140            indx = EXPphase[ODkey].split()
1141            ODkey = ODkey[:-1]+'B'
1142            vals = EXPphase[ODkey].split()
1143            for i,val in enumerate(vals):
1144                key = 'C(%s,%s,%s)'%(indx[3*i],indx[3*i+1],indx[3*i+2])
1145                shCoef[key] = float(val)
1146        textureData['SH Coeff'] = [False,shCoef]
1147       
1148    Phase = SetNewPhase(Name=PhaseName,SGData=SGData,cell=abc+angles+[Volume,])
1149    general = Phase['General']
1150    general['Type'] = Ptype
1151    general['SH Texture'] = textureData
1152    Phase['Atoms'] = Atoms
1153    return Phase
1154       
1155def ReadPDBPhase(filename):
1156    '''Read a phase from a PDB file.
1157    Called in the PDB phase import routine (see imports/G2phase.py)
1158    '''
1159    EightPiSq = 8.*math.pi**2
1160    file = open(filename, 'Ur')
1161    Phase = {}
1162    Title = ''
1163    Compnd = ''
1164    Atoms = []
1165    A = np.zeros(shape=(3,3))
1166    S = file.readline()
1167    while S:
1168        Atom = []
1169        if 'TITLE' in S[:5]:
1170            Title = S[10:72].strip()
1171            S = file.readline()
1172        elif 'COMPND    ' in S[:10]:
1173            Compnd = S[10:72].strip()
1174            S = file.readline()
1175        elif 'CRYST' in S[:5]:
1176            abc = S[7:34].split()
1177            angles = S[34:55].split()
1178            cell=[float(abc[0]),float(abc[1]),float(abc[2]),
1179                float(angles[0]),float(angles[1]),float(angles[2])]
1180            Volume = G2lat.calc_V(G2lat.cell2A(cell))
1181            AA,AB = G2lat.cell2AB(cell)
1182            SpGrp = S[55:65]
1183            E,SGData = G2spc.SpcGroup(SpGrp)
1184            # space group processing failed, try to look up name in table
1185            if E:
1186                SpGrpNorm = G2spc.StandardizeSpcName(SpGrp)
1187                if SpGrpNorm:
1188                    E,SGData = G2spc.SpcGroup(SpGrpNorm)
1189            while E:
1190                print G2spc.SGErrors(E)
1191                dlg = wx.TextEntryDialog(None,
1192                    SpGrp[:-1]+' is invalid \nN.B.: make sure spaces separate axial fields in symbol',
1193                    'ERROR in space group symbol','',style=wx.OK)
1194                if dlg.ShowModal() == wx.ID_OK:
1195                    SpGrp = dlg.GetValue()
1196                    E,SGData = G2spc.SpcGroup(SpGrp)
1197                else:
1198                    return None
1199                dlg.Destroy()               
1200            SGlines = G2spc.SGPrint(SGData)
1201            for line in SGlines: print line
1202            S = file.readline()
1203        elif 'SCALE' in S[:5]:
1204            V = (S[10:41].split())
1205            A[int(S[5])-1] = [float(V[0]),float(V[1]),float(V[2])]
1206            S = file.readline()
1207        elif 'ATOM' in S[:4] or 'HETATM' in S[:6]:
1208            XYZ = [float(S[31:39]),float(S[39:47]),float(S[47:55])]
1209            XYZ = np.inner(AB,XYZ)
1210            XYZ = np.where(abs(XYZ)<0.00001,0,XYZ)
1211            SytSym,Mult = G2spc.SytSym(XYZ,SGData)
1212            Uiso = float(S[61:67])/EightPiSq
1213            Type = S[12:14].lower()
1214            if Type[0] in '123456789':
1215                Type = Type[1:]
1216            Atom = [S[22:27].strip(),S[17:20].upper(),S[20:22],
1217                S[12:17].strip(),Type.strip().capitalize(),'',XYZ[0],XYZ[1],XYZ[2],
1218                float(S[55:61]),SytSym,Mult,'I',Uiso,0,0,0,0,0,0]
1219            S = file.readline()
1220            if 'ANISOU' in S[:6]:
1221                Uij = S[30:72].split()
1222                Uij = [float(Uij[0])/10000.,float(Uij[1])/10000.,float(Uij[2])/10000.,
1223                    float(Uij[3])/10000.,float(Uij[4])/10000.,float(Uij[5])/10000.]
1224                Atom = Atom[:14]+Uij
1225                Atom[12] = 'A'
1226                S = file.readline()
1227            Atom.append(ran.randint(0,sys.maxint))
1228            Atoms.append(Atom)
1229        else:           
1230            S = file.readline()
1231    file.close()
1232    if Title:
1233        PhaseName = Title
1234    elif Compnd:
1235        PhaseName = Compnd
1236    else:
1237        PhaseName = 'None'
1238    Phase = SetNewPhase(Name=PhaseName,SGData=SGData,cell=cell+[Volume,])
1239    Phase['General']['Type'] = 'macromolecular'
1240    Phase['Atoms'] = Atoms
1241   
1242    return Phase
1243
1244class MultipleChoicesDialog(wx.Dialog):
1245    '''A dialog that offers a series of choices, each with a
1246    title and a wx.Choice widget. Intended to be used Modally.
1247    typical input:
1248
1249        *  choicelist=[ ('a','b','c'), ('test1','test2'),('no choice',)]
1250        *  headinglist = [ 'select a, b or c', 'select 1 of 2', 'No option here']
1251       
1252    selections are placed in self.chosen when OK is pressed
1253    '''
1254    def __init__(self,choicelist,headinglist,
1255                 head='Select options',
1256                 title='Please select from options below',
1257                 parent=None):
1258        self.chosen = []
1259        wx.Dialog.__init__(
1260            self,parent,wx.ID_ANY,head, 
1261            pos=wx.DefaultPosition,style=wx.DEFAULT_DIALOG_STYLE)
1262        panel = wx.Panel(self)
1263        mainSizer = wx.BoxSizer(wx.VERTICAL)
1264        mainSizer.Add((10,10),1)
1265        topLabl = wx.StaticText(panel,wx.ID_ANY,title)
1266        mainSizer.Add(topLabl,0,wx.ALIGN_CENTER_VERTICAL|wx.CENTER,10)
1267        self.ChItems = []
1268        for choice,lbl in zip(choicelist,headinglist):
1269            mainSizer.Add((10,10),1)
1270            self.chosen.append(0)
1271            topLabl = wx.StaticText(panel,wx.ID_ANY,' '+lbl)
1272            mainSizer.Add(topLabl,0,wx.ALIGN_LEFT,10)
1273            self.ChItems.append(wx.Choice(self, wx.ID_ANY, (100, 50), choices = choice))
1274            mainSizer.Add(self.ChItems[-1],0,wx.ALIGN_CENTER,10)
1275
1276        OkBtn = wx.Button(panel,-1,"Ok")
1277        OkBtn.Bind(wx.EVT_BUTTON, self.OnOk)
1278        cancelBtn = wx.Button(panel,-1,"Cancel")
1279        cancelBtn.Bind(wx.EVT_BUTTON, self.OnCancel)
1280        btnSizer = wx.BoxSizer(wx.HORIZONTAL)
1281        btnSizer.Add((20,20),1)
1282        btnSizer.Add(OkBtn)
1283        btnSizer.Add((20,20),1)
1284        btnSizer.Add(cancelBtn)
1285        btnSizer.Add((20,20),1)
1286        mainSizer.Add(btnSizer,0,wx.EXPAND|wx.BOTTOM|wx.TOP, 10)
1287        panel.SetSizer(mainSizer)
1288        panel.Fit()
1289        self.Fit()
1290       
1291    def OnOk(self,event):
1292        parent = self.GetParent()
1293        if parent is not None: parent.Raise()
1294        # save the results from the choice widgets
1295        self.chosen = []
1296        for w in self.ChItems:
1297            self.chosen.append(w.GetSelection())
1298        self.EndModal(wx.ID_OK)             
1299           
1300    def OnCancel(self,event):
1301        parent = self.GetParent()
1302        if parent is not None: parent.Raise()
1303        self.chosen = []
1304        self.EndModal(wx.ID_CANCEL)             
1305           
1306def ExtractFileFromZip(filename, selection=None, confirmread=True,
1307                       confirmoverwrite=True, parent=None,
1308                       multipleselect=False):
1309    '''If the filename is a zip file, extract a file from that
1310    archive.
1311
1312    :param list Selection: used to predefine the name of the file
1313      to be extracted. Filename case and zip directory name are
1314      ignored in selection; the first matching file is used.
1315
1316    :param bool confirmread: if True asks the user to confirm before expanding
1317      the only file in a zip
1318
1319    :param bool confirmoverwrite: if True asks the user to confirm
1320      before overwriting if the extracted file already exists
1321
1322    :param bool multipleselect: if True allows more than one zip
1323      file to be extracted, a list of file(s) is returned.
1324      If only one file is present, do not ask which one, otherwise
1325      offer a list of choices (unless selection is used).
1326   
1327    :returns: the name of the file that has been created or a
1328      list of files (see multipleselect)
1329
1330    If the file is not a zipfile, return the name of the input file.
1331    If the zipfile is empty or no file has been selected, return None
1332    '''
1333    import zipfile # do this now, since we can save startup time by doing this only on need
1334    import shutil
1335    zloc = os.path.split(filename)[0]
1336    if not zipfile.is_zipfile(filename):
1337        #print("not zip")
1338        return filename
1339
1340    z = zipfile.ZipFile(filename,'r')
1341    zinfo = z.infolist()
1342
1343    if len(zinfo) == 0:
1344        #print('Zip has no files!')
1345        zlist = [-1]
1346    if selection:
1347        choices = [os.path.split(i.filename)[1].lower() for i in zinfo]
1348        if selection.lower() in choices:
1349            zlist = [choices.index(selection.lower())]
1350        else:
1351            print('debug: file '+str(selection)+' was not found in '+str(filename))
1352            zlist = [-1]
1353    elif len(zinfo) == 1 and confirmread:
1354        result = wx.ID_NO
1355        dlg = wx.MessageDialog(
1356            parent,
1357            'Is file '+str(zinfo[0].filename)+
1358            ' what you want to extract from '+
1359            str(os.path.split(filename)[1])+'?',
1360            'Confirm file', 
1361            wx.YES_NO | wx.ICON_QUESTION)
1362        try:
1363            result = dlg.ShowModal()
1364        finally:
1365            dlg.Destroy()
1366        if result == wx.ID_NO:
1367            zlist = [-1]
1368        else:
1369            zlist = [0]
1370    elif len(zinfo) == 1:
1371        zlist = [0]
1372    elif multipleselect:
1373        # select one or more from a from list
1374        choices = [i.filename for i in zinfo]
1375        dlg = wx.MultiChoiceDialog(parent,'Select file(s) to extract from zip file'+str(filename),
1376            'Choose file(s)',choices,wx.CHOICEDLG_STYLE,)
1377        if dlg.ShowModal() == wx.ID_OK:
1378            zlist = dlg.GetSelections()
1379        else:
1380            zlist = []
1381        dlg.Destroy()
1382    else:
1383        # select one from a from list
1384        choices = [i.filename for i in zinfo]
1385        dlg = wx.SingleChoiceDialog(parent,
1386            'Select file to extract from zip file'+str(filename),'Choose file',
1387            choices,)
1388        if dlg.ShowModal() == wx.ID_OK:
1389            zlist = [dlg.GetSelection()]
1390        else:
1391            zlist = [-1]
1392        dlg.Destroy()
1393       
1394    outlist = []
1395    for zindex in zlist:
1396        if zindex >= 0:
1397            efil = os.path.join(zloc, os.path.split(zinfo[zindex].filename)[1])
1398            if os.path.exists(efil) and confirmoverwrite:
1399                result = wx.ID_NO
1400                dlg = wx.MessageDialog(parent,
1401                    'File '+str(efil)+' already exists. OK to overwrite it?',
1402                    'Confirm overwrite',wx.YES_NO | wx.ICON_QUESTION)
1403                try:
1404                    result = dlg.ShowModal()
1405                finally:
1406                    dlg.Destroy()
1407                if result == wx.ID_NO:
1408                    zindex = -1
1409        if zindex >= 0:
1410            # extract the file to the current directory, regardless of it's original path
1411            #z.extract(zinfo[zindex],zloc)
1412            eloc,efil = os.path.split(zinfo[zindex].filename)
1413            outfile = os.path.join(zloc, efil)
1414            fpin = z.open(zinfo[zindex])
1415            fpout = file(outfile, "wb")
1416            shutil.copyfileobj(fpin, fpout)
1417            fpin.close()
1418            fpout.close()
1419            outlist.append(outfile)
1420    z.close()
1421    if multipleselect and len(outlist) >= 1:
1422        return outlist
1423    elif len(outlist) == 1:
1424        return outlist[0]
1425    else:
1426        return None
1427
1428######################################################################
1429# base classes for reading various types of data files
1430#   not used directly, only by subclassing
1431######################################################################
1432E,SGData = G2spc.SpcGroup('P 1') # data structure for default space group
1433class ImportBaseclass(object):
1434    '''Defines a base class for the reading of input files (diffraction
1435    data, coordinates,...
1436    '''
1437    def __init__(self,
1438                 formatName,
1439                 longFormatName=None,
1440                 extensionlist=[],
1441                 strictExtension=False,
1442                 ):
1443        self.formatName = formatName # short string naming file type
1444        if longFormatName: # longer string naming file type
1445            self.longFormatName = longFormatName
1446        else:
1447            self.longFormatName = formatName
1448        # define extensions that are allowed for the file type
1449        # for windows, remove any extensions that are duplicate, as case is ignored
1450        if sys.platform == 'windows' and extensionlist:
1451            extensionlist = list(set([s.lower() for s in extensionlist]))
1452        self.extensionlist = extensionlist
1453        # If strictExtension is True, the file will not be read, unless
1454        # the extension matches one in the extensionlist
1455        self.strictExtension = strictExtension
1456        self.warnings = ''
1457        # used for readers that will use multiple passes to read
1458        # more than one data block
1459        self.repeat = False
1460        self.repeatcount = 0
1461        #print 'created',self.__class__
1462
1463    def BlockSelector(self, ChoiceList, ParentFrame=None,
1464                      title='Select a block',
1465                      size=None, header='Block Selector',
1466                      useCancel=True):
1467        ''' Provide a wx dialog to select a block if the file contains more
1468        than one set of data and one must be selected
1469        '''
1470        if useCancel:
1471            dlg = wx.SingleChoiceDialog(
1472                ParentFrame,title, header,ChoiceList)
1473        else:
1474            dlg = wx.SingleChoiceDialog(
1475                ParentFrame,title, header,ChoiceList,
1476                style=wx.DEFAULT_DIALOG_STYLE|wx.RESIZE_BORDER|wx.OK|wx.CENTRE)
1477        if size: dlg.SetSize(size)
1478        dlg.CenterOnParent()
1479        if dlg.ShowModal() == wx.ID_OK:
1480            sel = dlg.GetSelection()
1481            return sel
1482        else:
1483            return None
1484        dlg.Destroy()
1485
1486    def MultipleBlockSelector(self, ChoiceList, ParentFrame=None,
1487        title='Select a block',size=None, header='Block Selector'):
1488        '''Provide a wx dialog to select a block of data if the
1489        file contains more than one set of data and one must be
1490        selected.
1491
1492        :returns: a list of the selected blocks
1493        '''
1494        dlg = wx.MultiChoiceDialog(ParentFrame,title, header,ChoiceList+['Select all'],
1495            wx.CHOICEDLG_STYLE)
1496        dlg.CenterOnParent()
1497        if size: dlg.SetSize(size)
1498        if dlg.ShowModal() == wx.ID_OK:
1499            sel = dlg.GetSelections()
1500        else:
1501            return []
1502        dlg.Destroy()
1503        selected = []
1504        if len(ChoiceList) in sel:
1505            return range(len(ChoiceList))
1506        else:
1507            return sel
1508        return selected
1509
1510    def MultipleChoicesDialog(self, choicelist, headinglist, ParentFrame=None, **kwargs):
1511        '''A modal dialog that offers a series of choices, each with a title and a wx.Choice
1512        widget. Typical input:
1513       
1514           * choicelist=[ ('a','b','c'), ('test1','test2'),('no choice',)]
1515           
1516           * headinglist = [ 'select a, b or c', 'select 1 of 2', 'No option here']
1517           
1518        optional keyword parameters are: head (window title) and title
1519        returns a list of selected indicies for each choice (or None)
1520        '''
1521        result = None
1522        dlg = MultipleChoicesDialog(choicelist,headinglist,
1523            parent=ParentFrame, **kwargs)         
1524        dlg.CenterOnParent()
1525        if dlg.ShowModal() == wx.ID_OK:
1526            result = dlg.chosen
1527        dlg.Destroy()
1528        return result
1529
1530    def ShowBusy(self):
1531        wx.BeginBusyCursor()
1532
1533    def DoneBusy(self):
1534        wx.EndBusyCursor()
1535       
1536#    def Reader(self, filename, filepointer, ParentFrame=None, **unused):
1537#        '''This method must be supplied in the child class
1538#        it will read the file
1539#        '''
1540#        return True # if read OK
1541#        return False # if an error occurs
1542
1543    def ExtensionValidator(self, filename):
1544        '''This methods checks if the file has the correct extension
1545        Return False if this filename will not be supported by this reader
1546        Return True if the extension matches the list supplied by the reader
1547        Return None if the reader allows un-registered extensions
1548        '''
1549        if filename:
1550            ext = os.path.splitext(filename)[1]
1551            if sys.platform == 'windows': ext = ext.lower()
1552            if ext in self.extensionlist: return True
1553            if self.strictExtension: return False
1554        return None
1555
1556    def ContentsValidator(self, filepointer):
1557        '''This routine will attempt to determine if the file can be read
1558        with the current format.
1559        This will typically be overridden with a method that
1560        takes a quick scan of [some of]
1561        the file contents to do a "sanity" check if the file
1562        appears to match the selected format.
1563        Expected to be called via self.Validator()
1564        '''
1565        #filepointer.seek(0) # rewind the file pointer
1566        return True
1567
1568class ImportPhase(ImportBaseclass):
1569    '''Defines a base class for the reading of files with coordinates
1570    '''
1571    def __init__(self,formatName,longFormatName=None,extensionlist=[],
1572        strictExtension=False,):
1573        # call parent __init__
1574        ImportBaseclass.__init__(self,formatName,longFormatName,
1575            extensionlist,strictExtension)
1576        # define a default Phase structure
1577        self.Phase = SetNewPhase(Name='new phase',SGData=SGData)
1578
1579    def PhaseSelector(self, ChoiceList, ParentFrame=None,
1580        title='Select a phase', size=None,header='Phase Selector'):
1581        ''' Provide a wx dialog to select a phase if the file contains more
1582        than one phase
1583        '''
1584        return self.BlockSelector(ChoiceList,ParentFrame,title,
1585            size,header)
1586
1587class ImportStructFactor(ImportBaseclass):
1588    '''Defines a base class for the reading of files with tables
1589    of structure factors
1590
1591    Note that the default controls are stored in self.Controls and the
1592    default instrument parameters are stored in self.Parameters.
1593    These can be changed, but any changes will be the defaults for all
1594    subsequent uses of the :class:`ImportStructFactor` derived classes
1595    until :meth:`InitControls` and :meth:`InitParameters` are
1596    called. Probably better to use :meth:`UpdateControls` and
1597    :meth:`UpdateControls` (adding new args if needed) to change
1598    values.
1599    '''
1600    def __init__(self,formatName,longFormatName=None,extensionlist=[],
1601        strictExtension=False,):
1602        ImportBaseclass.__init__(self,formatName,longFormatName,
1603            extensionlist,strictExtension)
1604
1605        # define contents of Structure Factor entry
1606        self.InitParameters()
1607        self.InitControls()
1608        self.RefDict = {'RefList':[],'FF':[]}
1609       
1610    def InitControls(self):
1611        'initialize the controls structure'
1612        self.Controls = { # dictionary with plotting controls
1613            'Type' : 'Fosq',
1614            'ifFc' : False,    #
1615            'HKLmax' : [None,None,None],
1616            'HKLmin' : [None,None,None],
1617            'FoMax' : None,   # maximum observed structure factor as Fo
1618            'Zone' : '001',
1619            'Layer' : 0,
1620            'Scale' : 1.0,
1621            'log-lin' : 'lin',
1622            }
1623
1624    def InitParameters(self):
1625        'initialize the instrument parameters structure'
1626        Lambda = 0.70926
1627        HistType = 'SXC'
1628        self.Parameters = [{'Type':[HistType,HistType], # create the structure
1629                            'Lam':[Lambda,Lambda]
1630                            }, {}]
1631
1632    def UpdateParameters(self,Type=None,Wave=None):
1633        'Revise the instrument parameters'
1634        if Type is not None:
1635            self.Parameters[0]['Type'] = [Type,Type]
1636        if Wave is not None:
1637            self.Parameters[0]['Lam'] = [Wave,Wave]
1638           
1639    def UpdateControls(self,Type='Fosq',FcalcPresent=False):
1640        '''Scan through the reflections to update the Controls dictionary
1641        '''
1642        self.Controls['Type'] = Type
1643        self.Controls['ifFc'] = FcalcPresent
1644        HKLmax = [None,None,None]
1645        HKLmin = [None,None,None]
1646        Fo2max = None
1647        for refl in self.RefDict['RefList']:
1648            HKL = refl[:3]
1649            if Fo2max is None:
1650                Fo2max = refl[8]
1651            else:
1652                Fo2max = max(Fo2max,refl[8])
1653            for i,hkl in enumerate(HKL):
1654                if HKLmax[i] is None:
1655                    HKLmax[i] = hkl
1656                    HKLmin[i] = hkl
1657                else:
1658                    HKLmax[i] = max(HKLmax[i],hkl)
1659                    HKLmin[i] = min(HKLmin[i],hkl)
1660        self.Controls['HKLmax'] = HKLmax
1661        self.Controls['HKLmin'] = HKLmin
1662        if Type ==  'Fosq':
1663            self.Controls['FoMax'] = np.sqrt(Fo2max)
1664        elif Type ==  'Fo':
1665            self.Controls['FoMax'] = Fo2max
1666        else:
1667            print "Unsupported Struct Fact type in ImportStructFactor.UpdateControls"
1668            raise Exception,"Unsupported Struct Fact type in ImportStructFactor.UpdateControls"
1669
1670######################################################################
1671class ImportPowderData(ImportBaseclass):
1672    '''Defines a base class for the reading of files with powder data
1673    '''
1674    # define some default instrument parameter files
1675    # just like GSAS, sigh
1676    defaultIparm_lbl = []
1677    defaultIparms = []
1678    defaultIparm_lbl.append('CuKa lab data')
1679    defaultIparms.append({
1680        'INS   HTYPE ':'PXC ',
1681        'INS  1 ICONS':'  1.540500  1.544300       0.0         0       0.7    0       0.5   ',
1682        'INS  1PRCF1 ':'    3    8      0.01                                                ',
1683        'INS  1PRCF11':'   2.000000E+00  -2.000000E+00   5.000000E+00   0.000000E+00        ',
1684        'INS  1PRCF12':'   0.000000E+00   0.000000E+00   0.150000E-01   0.150000E-01        ',
1685        })
1686    defaultIparm_lbl.append('0.6A synch')
1687    defaultIparms.append({
1688        'INS   HTYPE ':'PXC ',
1689        'INS  1 ICONS':'  0.600000  0.000000       0.0         0      0.99    0       0.5   ',
1690        'INS  1PRCF1 ':'    3    8      0.01                                                ',
1691        'INS  1PRCF11':'   1.000000E+00  -1.000000E+00   0.300000E+00   0.000000E+00        ',
1692        'INS  1PRCF12':'   0.000000E+00   0.000000E+00   0.100000E-01   0.100000E-01        ',
1693        })
1694    defaultIparm_lbl.append('1.5A CW neutron data')
1695    defaultIparms.append({
1696        'INS   HTYPE ':'PNC',
1697        'INS  1 ICONS':'   1.54020   0.00000   0.04000         0',
1698        'INS  1PRCF1 ':'    3    8     0.005',
1699        'INS  1PRCF11':'   0.239700E+03  -0.298200E+03   0.180800E+03   0.000000E+00',
1700        'INS  1PRCF12':'   0.000000E+00   0.000000E+00   0.400000E-01   0.300000E-01',
1701        })
1702    defaultIparm_lbl.append('10m TOF backscattering bank')
1703    defaultIparms.append({
1704        'INS   HTYPE ':'PNT',
1705        'INS  1 ICONS':'   5000.00      0.00      0.00',
1706        'INS  1BNKPAR':'    1.0000   150.000',       
1707        'INS  1PRCF1 ':'    1    8   0.01000',
1708        'INS  1PRCF11':'   0.000000E+00   5.000000E+00   3.000000E-02   1.000000E-03',
1709        'INS  1PRCF12':'   0.000000E+00   4.000000E+01   0.000000E+00   0.000000E+00',       
1710        })
1711    defaultIparm_lbl.append('10m TOF 90deg bank')
1712    defaultIparms.append({
1713        'INS   HTYPE ':'PNT',
1714        'INS  1 ICONS':'   3500.00      0.00      0.00',
1715        'INS  1BNKPAR':'    1.0000    90.000',       
1716        'INS  1PRCF1 ':'    1    8   0.01000',
1717        'INS  1PRCF11':'   0.000000E+00   5.000000E+00   3.000000E-02   4.000000E-03',
1718        'INS  1PRCF12':'   0.000000E+00   8.000000E+01   0.000000E+00   0.000000E+00',       
1719        })
1720    defaultIparm_lbl.append('63m POWGEN 90deg bank')
1721    defaultIparms.append({
1722        'INS   HTYPE ':'PNT',
1723        'INS  1 ICONS':'  22585.80      0.00      0.00',
1724        'INS  1BNKPAR':'    1.0000    90.000',       
1725        'INS  1PRCF1 ':'    1    8   0.01000',
1726        'INS  1PRCF11':'   0.000000E+00   1.000000E+00   3.000000E-02   4.000000E-03',
1727        'INS  1PRCF12':'   0.000000E+00   8.000000E+01   0.000000E+00   0.000000E+00',       
1728        })
1729    def __init__(self,
1730                 formatName,
1731                 longFormatName=None,
1732                 extensionlist=[],
1733                 strictExtension=False,
1734                 ):
1735        ImportBaseclass.__init__(self,formatName,
1736                                            longFormatName,
1737                                            extensionlist,
1738                                            strictExtension)
1739        self.powderentry = ['',None,None] #  (filename,Pos,Bank)
1740        self.powderdata = [] # Powder dataset
1741        '''A powder data set is a list with items [x,y,w,yc,yb,yd]:
1742                np.array(x), # x-axis values
1743                np.array(y), # powder pattern intensities
1744                np.array(w), # 1/sig(intensity)^2 values (weights)
1745                np.array(yc), # calc. intensities (zero)
1746                np.array(yb), # calc. background (zero)
1747                np.array(yd), # obs-calc profiles
1748        '''                           
1749        self.comments = []
1750        self.idstring = ''
1751        self.Sample = G2pdG.SetDefaultSample()
1752        self.GSAS = None     # used in TOF
1753        self.clockWd = None  # used in TOF
1754        self.repeat_instparm = True # Should a parm file be
1755        #                             used for multiple histograms?
1756        self.instparm = None # name hint
1757        self.instfile = '' # full path name to instrument parameter file
1758        self.instbank = '' # inst parm bank number
1759        self.instmsg = ''  # a label that gets printed to show
1760                           # where instrument parameters are from
1761        self.numbanks = 1
1762        self.instdict = {} # place items here that will be transferred to the instrument parameters
1763######################################################################
1764class ExportBaseclass(object):
1765    '''Defines a base class for the exporting of GSAS-II results
1766    '''
1767    def __init__(self,
1768                 G2frame,
1769                 formatName,
1770                 extension,
1771                 longFormatName=None,
1772                 ):
1773        self.G2frame = G2frame
1774        self.formatName = formatName # short string naming file type
1775        self.extension = extension
1776        if longFormatName: # longer string naming file type
1777            self.longFormatName = longFormatName
1778        else:
1779            self.longFormatName = formatName
1780        self.OverallParms = {}
1781        self.Phases = {}
1782        self.Histograms = {}
1783        self.powderDict = {}
1784        self.xtalDict = {}
1785        self.parmDict = {}
1786        self.sigDict = {}
1787        # updated in InitExport:
1788        self.currentExportType = None # type of export that has been requested
1789        # updated in ExportSelect (when used):
1790        self.phasenam = None # a list of selected phases
1791        self.histnam = None # a list of selected histograms
1792        self.filename = None # name of file to be written
1793       
1794        # items that should be defined in a subclass of this class
1795        self.exporttype = []  # defines the type(s) of exports that the class can handle.
1796        # The following types are defined: 'project', "phase", "powder", "single"
1797        self.multiple = False # set as True if the class can export multiple phases or histograms
1798        # self.multiple is ignored for "project" exports
1799
1800    def InitExport(self,event):
1801        '''Determines the type of menu that called the Exporter.
1802        '''
1803        if event:
1804            self.currentExportType = self.G2frame.ExportLookup.get(event.Id)
1805
1806    def ExportSelect(self,AskFile=True):
1807        '''Selects histograms or phases when needed. Sets a default file name.
1808
1809        :param bool AskFile: if AskFile is True (default) get the name of the file
1810          in a dialog
1811        :returns: True in case of an error
1812        '''
1813       
1814        if self.currentExportType == 'phase':
1815            if len(self.Phases) == 0:
1816                self.G2frame.ErrorDialog(
1817                    'Empty project',
1818                    'Project does not contain any phases.')
1819                return True
1820            elif len(self.Phases) == 1:
1821                self.phasenam = self.Phases.keys()
1822            elif self.multiple: 
1823                choices = sorted(self.Phases.keys())
1824                phasenum = G2gd.ItemSelector(choices,self.G2frame,multiple=True)
1825                if phasenum is None: return True
1826                self.phasenam = [choices[i] for i in phasenum]
1827                if not self.phasenam: return True
1828            else:
1829                choices = sorted(self.Phases.keys())
1830                phasenum = G2gd.ItemSelector(choices,self.G2frame)
1831                if phasenum is None: return True
1832                self.phasenam = [choices[phasenum]]
1833        elif self.currentExportType == 'single':
1834            if len(self.xtalDict) == 0:
1835                self.G2frame.ErrorDialog(
1836                    'Empty project',
1837                    'Project does not contain any single crystal data.')
1838                return True
1839            elif len(self.xtalDict) == 1:
1840                self.histnam = self.xtalDict.values()
1841            elif self.multiple:
1842                choices = sorted(self.xtalDict.values())
1843                hnum = G2gd.ItemSelector(choices,self.G2frame,multiple=True)
1844                if not hnum: return True
1845                self.histnam = [choices[i] for i in hnum]
1846            else:
1847                choices = sorted(self.xtalDict.values())
1848                hnum = G2gd.ItemSelector(choices,self.G2frame)
1849                if hnum is None: return True
1850                self.histnam = [choices[hnum]]
1851        elif self.currentExportType == 'powder':
1852            if len(self.powderDict) == 0:
1853                self.G2frame.ErrorDialog(
1854                    'Empty project',
1855                    'Project does not contain any powder data.')
1856                return True
1857            elif len(self.powderDict) == 1:
1858                self.histnam = self.powderDict.values()
1859            elif self.multiple:
1860                choices = sorted(self.powderDict.values())
1861                hnum = G2gd.ItemSelector(choices,self.G2frame,multiple=True)
1862                if not hnum: return True
1863                self.histnam = [choices[i] for i in hnum]
1864            else:
1865                choices = sorted(self.powderDict.values())
1866                hnum = G2gd.ItemSelector(choices,self.G2frame)
1867                if hnum is None: return True
1868                self.histnam = [choices[hnum]]
1869        elif self.currentExportType == 'image':
1870            if len(self.Histograms) == 0:
1871                self.G2frame.ErrorDialog(
1872                    'Empty project',
1873                    'Project does not contain any images.')
1874                return True
1875            elif len(self.Histograms) == 1:
1876                self.histnam = self.Histograms.keys()
1877            else:
1878                choices = sorted(self.Histograms.keys())
1879                hnum = G2gd.ItemSelector(choices,self.G2frame,multiple=self.multiple)
1880                if self.multiple:
1881                    if not hnum: return True
1882                    self.histnam = [choices[i] for i in hnum]
1883                else:
1884                    if hnum is None: return True
1885                    self.histnam = [choices[hnum]]
1886        if self.currentExportType == 'map':
1887            # search for phases with maps
1888            mapPhases = []
1889            choices = []
1890            for phasenam in sorted(self.Phases):
1891                phasedict = self.Phases[phasenam] # pointer to current phase info           
1892                if len(phasedict['General']['Map'].get('rho',[])):
1893                    mapPhases.append(phasenam)
1894                    if phasedict['General']['Map'].get('Flip'):
1895                        choices.append('Charge flip map: '+str(phasenam))
1896                    elif phasedict['General']['Map'].get('MapType'):
1897                        choices.append(
1898                            str(phasedict['General']['Map'].get('MapType'))
1899                            + ' map: ' + str(phasenam))
1900                    else:
1901                        choices.append('unknown map: '+str(phasenam))
1902            # select a map if needed
1903            if len(mapPhases) == 0:
1904                self.G2frame.ErrorDialog(
1905                    'Empty project',
1906                    'Project does not contain any maps.')
1907                return True
1908            elif len(mapPhases) == 1:
1909                self.phasenam = mapPhases
1910            else: 
1911                phasenum = G2gd.ItemSelector(choices,self.G2frame,multiple=self.multiple)
1912                if self.multiple:
1913                    if not phasenum: return True
1914                    self.phasenam = [mapPhases[i] for i in phasenum]
1915                else:
1916                    if phasenum is None: return True
1917                    self.phasenam = [mapPhases[phasenum]]
1918
1919        if AskFile:
1920            self.filename = self.askSaveFile()
1921        else:
1922            self.filename = self.defaultSaveFile()
1923        if not self.filename: return True
1924       
1925    # def SetupExport(self,event,AskFile=True):
1926    #     '''Determines the type of menu that called the Exporter. Selects histograms
1927    #     or phases when needed. Better to replace with individual calls to
1928    #     self.InitExport() and self.ExportSelect() so that the former can be called prior
1929    #     to self.LoadTree()
1930
1931    #     :param bool AskFile: if AskFile is True (default) get the name of the file
1932    #       in a dialog
1933    #     :returns: True in case of an error
1934    #     '''
1935    #     self.ExportInit(event)
1936    #     return self.ExportSelect(AskFile)
1937                   
1938    def loadParmDict(self):
1939        '''Load the GSAS-II refinable parameters from the tree into a dict (self.parmDict). Update
1940        refined values to those from the last cycle and set the uncertainties for the
1941        refined parameters in another dict (self.sigDict).
1942
1943        Expands the parm & sig dicts to include values derived from constraints.
1944        '''
1945        self.parmDict = {}
1946        self.sigDict = {}
1947        rigidbodyDict = {}
1948        covDict = {}
1949        consDict = {}
1950        Histograms,Phases = self.G2frame.GetUsedHistogramsAndPhasesfromTree()
1951        if self.G2frame.PatternTree.IsEmpty(): return # nothing to do
1952        item, cookie = self.G2frame.PatternTree.GetFirstChild(self.G2frame.root)
1953        while item:
1954            name = self.G2frame.PatternTree.GetItemText(item)
1955            if name == 'Rigid bodies':
1956                 rigidbodyDict = self.G2frame.PatternTree.GetItemPyData(item)
1957            elif name == 'Covariance':
1958                 covDict = self.G2frame.PatternTree.GetItemPyData(item)
1959            elif name == 'Constraints':
1960                 consDict = self.G2frame.PatternTree.GetItemPyData(item)
1961            item, cookie = self.G2frame.PatternTree.GetNextChild(self.G2frame.root, cookie)
1962        rbVary,rbDict =  G2stIO.GetRigidBodyModels(rigidbodyDict,Print=False)
1963        self.parmDict.update(rbDict)
1964        rbIds = rigidbodyDict.get('RBIds',{'Vector':[],'Residue':[]})
1965        Natoms,atomIndx,phaseVary,phaseDict,pawleyLookup,FFtables,BLtables =  G2stIO.GetPhaseData(
1966            Phases,RestraintDict=None,rbIds=rbIds,Print=False)
1967        self.parmDict.update(phaseDict)
1968        hapVary,hapDict,controlDict =  G2stIO.GetHistogramPhaseData(
1969            Phases,Histograms,Print=False,resetRefList=False)
1970        self.parmDict.update(hapDict)
1971        histVary,histDict,controlDict =  G2stIO.GetHistogramData(Histograms,Print=False)
1972        self.parmDict.update(histDict)
1973        self.parmDict.update(zip(
1974            covDict.get('varyList',[]),
1975            covDict.get('variables',[])))
1976        self.sigDict = dict(zip(
1977            covDict.get('varyList',[]),
1978            covDict.get('sig',[])))
1979        # expand to include constraints: first compile a list of constraints
1980        constList = []
1981        for item in consDict:
1982            constList += consDict[item]
1983        # now process the constraints
1984        G2mv.InitVars()
1985        constDict,fixedList,ignored = G2stIO.ProcessConstraints(constList)
1986        varyList = covDict.get('varyListStart')
1987        if varyList is None and len(constDict) == 0:
1988            # no constraints can use varyList
1989            varyList = covDict.get('varyList')
1990        elif varyList is None:
1991            # old GPX file from before pre-constraint varyList is saved
1992            print ' *** Old refinement: Please use Calculate/Refine to redo  ***'
1993            raise Exception(' *** Export aborted ***')
1994        else:
1995            varyList = list(varyList)
1996        try:
1997            groups,parmlist = G2mv.GroupConstraints(constDict)
1998            G2mv.GenerateConstraints(groups,parmlist,varyList,constDict,fixedList)
1999        except:
2000            # this really should not happen
2001            print ' *** ERROR - constraints are internally inconsistent ***'
2002            errmsg, warnmsg = G2mv.CheckConstraints(varyList,constDict,fixedList)
2003            print 'Errors',errmsg
2004            if warnmsg: print 'Warnings',warnmsg
2005            raise Exception(' *** CIF creation aborted ***')
2006        # add the constrained values to the parameter dictionary
2007        G2mv.Dict2Map(self.parmDict,varyList)
2008        # and add their uncertainties into the esd dictionary (sigDict)
2009        if covDict.get('covMatrix') is not None:
2010            self.sigDict.update(G2mv.ComputeDepESD(covDict['covMatrix'],covDict['varyList'],self.parmDict))
2011
2012    def loadTree(self):
2013        '''Load the contents of the data tree into a set of dicts
2014        (self.OverallParms, self.Phases and self.Histogram as well as self.powderDict
2015        & self.xtalDict)
2016       
2017        * The childrenless data tree items are overall parameters/controls for the
2018          entire project and are placed in self.OverallParms
2019        * Phase items are placed in self.Phases
2020        * Data items are placed in self.Histogram. The key for these data items
2021          begin with a keyword, such as PWDR, IMG, HKLF,... that identifies the data type.
2022        '''
2023        self.OverallParms = {}
2024        self.powderDict = {}
2025        self.xtalDict = {}
2026        if self.G2frame.PatternTree.IsEmpty(): return # nothing to do
2027        histType = None       
2028        if self.currentExportType == 'phase':
2029            # if exporting phases load them here
2030            sub = G2gd.GetPatternTreeItemId(self.G2frame,self.G2frame.root,'Phases')
2031            if not sub:
2032                print 'no phases found'
2033                return True
2034            item, cookie = self.G2frame.PatternTree.GetFirstChild(sub)
2035            while item:
2036                phaseName = self.G2frame.PatternTree.GetItemText(item)
2037                self.Phases[phaseName] =  self.G2frame.PatternTree.GetItemPyData(item)
2038                item, cookie = self.G2frame.PatternTree.GetNextChild(sub, cookie)
2039            return
2040        elif self.currentExportType == 'single':
2041            histType = 'HKLF'
2042        elif self.currentExportType == 'powder':
2043            histType = 'PWDR'
2044        elif self.currentExportType == 'image':
2045            histType = 'IMG'
2046
2047        if histType: # Loading just one kind of tree entry
2048            item, cookie = self.G2frame.PatternTree.GetFirstChild(self.G2frame.root)
2049            while item:
2050                name = self.G2frame.PatternTree.GetItemText(item)
2051                if name.startswith(histType):
2052                    if self.Histograms.get(name): # there is already an item with this name
2053                        if name[-1] == '9':
2054                            name = name[:-1] + '10'
2055                        elif name[-1] in '012345678':
2056                            name = name[:-1] + str(int(name[:-1])+1)
2057                        else:                           
2058                            name += '-1'
2059                    self.Histograms[name] = {}
2060                    # the main info goes into Data, but the 0th
2061                    # element contains refinement results, carry
2062                    # that over too now.
2063                    self.Histograms[name]['Data'] = self.G2frame.PatternTree.GetItemPyData(item)[1]
2064                    self.Histograms[name][0] = self.G2frame.PatternTree.GetItemPyData(item)[0]
2065                    item2, cookie2 = self.G2frame.PatternTree.GetFirstChild(item)
2066                    while item2: 
2067                        child = self.G2frame.PatternTree.GetItemText(item2)
2068                        self.Histograms[name][child] = self.G2frame.PatternTree.GetItemPyData(item2)
2069                        item2, cookie2 = self.G2frame.PatternTree.GetNextChild(item, cookie2)
2070                item, cookie = self.G2frame.PatternTree.GetNextChild(self.G2frame.root, cookie)
2071            # index powder and single crystal histograms by number
2072            for hist in self.Histograms:
2073                if hist.startswith("PWDR"): 
2074                    d = self.powderDict
2075                elif hist.startswith("HKLF"): 
2076                    d = self.xtalDict
2077                else:
2078                    return                   
2079                i = self.Histograms[hist].get('hId')
2080                if i is None and not d.keys():
2081                    i = 0
2082                elif i is None or i in d.keys():
2083                    i = max(d.keys())+1
2084                d[i] = hist
2085            return
2086        # else standard load: using all interlinked phases and histograms
2087        self.Histograms,self.Phases = self.G2frame.GetUsedHistogramsAndPhasesfromTree()
2088        item, cookie = self.G2frame.PatternTree.GetFirstChild(self.G2frame.root)
2089        while item:
2090            name = self.G2frame.PatternTree.GetItemText(item)
2091            item2, cookie2 = self.G2frame.PatternTree.GetFirstChild(item)
2092            if not item2: 
2093                self.OverallParms[name] = self.G2frame.PatternTree.GetItemPyData(item)
2094            item, cookie = self.G2frame.PatternTree.GetNextChild(self.G2frame.root, cookie)
2095        # index powder and single crystal histograms
2096        for hist in self.Histograms:
2097            i = self.Histograms[hist]['hId']
2098            if hist.startswith("PWDR"): 
2099                self.powderDict[i] = hist
2100            elif hist.startswith("HKLF"): 
2101                self.xtalDict[i] = hist
2102
2103    def dumpTree(self,mode='type'):
2104        '''Print out information on the data tree dicts loaded in loadTree
2105        '''
2106        print '\nOverall'
2107        if mode == 'type':
2108            def Show(arg): return type(arg)
2109        else:
2110            def Show(arg): return arg
2111        for key in self.OverallParms:
2112            print '  ',key,Show(self.OverallParms[key])
2113        print 'Phases'
2114        for key1 in self.Phases:
2115            print '    ',key1,Show(self.Phases[key1])
2116        print 'Histogram'
2117        for key1 in self.Histograms:
2118            print '    ',key1,Show(self.Histograms[key1])
2119            for key2 in self.Histograms[key1]:
2120                print '      ',key2,Show(self.Histograms[key1][key2])
2121
2122    def defaultSaveFile(self):
2123        return os.path.abspath(
2124            os.path.splitext(self.G2frame.GSASprojectfile
2125                             )[0]+self.extension)
2126       
2127    def askSaveFile(self):
2128        '''Ask the user to supply a file name
2129
2130        :returns: a file name (str)
2131        '''
2132        defnam = os.path.splitext(
2133            os.path.split(self.G2frame.GSASprojectfile)[1]
2134            )[0]+self.extension
2135        dlg = wx.FileDialog(
2136            self.G2frame, 'Input name for file to write', '.', defnam,
2137            self.longFormatName+' (*'+self.extension+')|*'+self.extension,
2138            wx.FD_SAVE|wx.FD_OVERWRITE_PROMPT|wx.CHANGE_DIR)
2139        dlg.CenterOnParent()
2140        try:
2141            if dlg.ShowModal() == wx.ID_OK:
2142                filename = dlg.GetPath()
2143                # make sure extension is correct
2144                filename = os.path.splitext(filename)[0]+self.extension
2145            else:
2146                filename = None
2147        finally:
2148            dlg.Destroy()
2149        return filename
2150
2151    # Tools for file writing.
2152    def OpenFile(self,fil=None):
2153        '''Open the output file
2154
2155        :param str fil: The name of the file to open. If None (default)
2156          the name defaults to self.filename.
2157        :returns: the file object opened by the routine which is also
2158          saved as self.fp
2159        '''
2160        if not fil:
2161            fil = self.filename
2162        self.fp = open(fil,'w')
2163        return self.fp
2164    def Write(self,line):
2165        '''write a line of output, attaching a line-end character
2166
2167        :param str line: the text to be written.
2168        '''
2169        self.fp.write(line+'\n')
2170    def CloseFile(self,fp=None):
2171        '''Close a file opened in OpenFile
2172
2173        :param file fp: the file object to be closed. If None (default)
2174          file object self.fp is closed.
2175        '''
2176        if fp is None:
2177            fp = self.fp
2178            self.fp = None
2179        fp.close()
2180    # Tools to pull information out of the data arrays
2181    def GetCell(self,phasenam):
2182        """Gets the unit cell parameters and their s.u.'s for a selected phase
2183
2184        :param str phasenam: the name for the selected phase
2185        :returns: `cellList,cellSig` where each is a 7 element list corresponding
2186          to a, b, c, alpha, beta, gamma, volume where `cellList` has the
2187          cell values and `cellSig` has their uncertainties.
2188        """
2189        phasedict = self.Phases[phasenam] # pointer to current phase info
2190        try:
2191            pfx = str(phasedict['pId'])+'::'
2192            A,sigA = G2stIO.cellFill(pfx,phasedict['General']['SGData'],self.parmDict,self.sigDict)
2193            cellSig = G2stIO.getCellEsd(pfx,
2194                                        phasedict['General']['SGData'],A,
2195                                        self.OverallParms['Covariance'])  # returns 7 vals, includes sigVol
2196            cellList = G2lat.A2cell(A) + (G2lat.calc_V(A),)
2197            return cellList,cellSig
2198        except KeyError:
2199            cell = phasedict['General']['Cell'][1:]
2200            return cell,7*[0]
2201   
2202    def GetAtoms(self,phasenam):
2203        """Gets the atoms associated with a phase. Can be used with standard
2204        or macromolecular phases
2205
2206        :param str phasenam: the name for the selected phase
2207        :returns: a list of items for eac atom where each item is a list containing:
2208          label, typ, mult, xyz, and td, where
2209
2210          * label and typ are the atom label and the scattering factor type (str)
2211          * mult is the site multiplicity (int)
2212          * xyz is contains a list with four pairs of numbers:
2213            x, y, z and fractional occupancy and
2214            their standard uncertainty (or a negative value)
2215          * td is contains a list with either one or six pairs of numbers:
2216            if one number it is U\ :sub:`iso` and with six numbers it is
2217            U\ :sub:`11`, U\ :sub:`22`, U\ :sub:`33`, U\ :sub:`12`, U\ :sub:`13` & U\ :sub:`23`
2218            paired with their standard uncertainty (or a negative value)
2219        """
2220        phasedict = self.Phases[phasenam] # pointer to current phase info           
2221        cx,ct,cs,cia = phasedict['General']['AtomPtrs']
2222        cfrac = cx+3
2223        fpfx = str(phasedict['pId'])+'::Afrac:'       
2224        atomslist = []
2225        for i,at in enumerate(phasedict['Atoms']):
2226            if phasedict['General']['Type'] == 'macromolecular':
2227                label = '%s_%s_%s_%s'%(at[ct-1],at[ct-3],at[ct-4],at[ct-2])
2228            else:
2229                label = at[ct-1]
2230            fval = self.parmDict.get(fpfx+str(i),at[cfrac])
2231            fsig = self.sigDict.get(fpfx+str(i),-0.009)
2232            mult = at[cs+1]
2233            typ = at[ct]
2234            xyz = []
2235            for j,v in enumerate(('x','y','z')):
2236                val = at[cx+j]
2237                pfx = str(phasedict['pId'])+'::dA'+v+':'+str(i)
2238                sig = self.sigDict.get(pfx,-0.000009)
2239                xyz.append((val,sig))
2240            xyz.append((fval,fsig))
2241            td = []
2242            if at[cia] == 'I':
2243                pfx = str(phasedict['pId'])+'::AUiso:'+str(i)
2244                val = self.parmDict.get(pfx,at[cia+1])
2245                sig = self.sigDict.get(pfx,-0.0009)
2246                td.append((val,sig))
2247            else:
2248                for i,var in enumerate(('AU11','AU22','AU33','AU12','AU13','AU23')):
2249                    pfx = str(phasedict['pId'])+'::'+var+':'+str(i)
2250                    val = self.parmDict.get(pfx,at[cia+2+i])
2251                    sig = self.sigDict.get(pfx,-0.0009)
2252                    td.append((val,sig))
2253            atomslist.append((label,typ,mult,xyz,td))
2254        return atomslist
2255######################################################################
2256
2257def ReadCIF(URLorFile):
2258    '''Open a CIF, which may be specified as a file name or as a URL using PyCifRW
2259    (from James Hester).
2260    The open routine gets confused with DOS names that begin with a letter and colon
2261    "C:\dir\" so this routine will try to open the passed name as a file and if that
2262    fails, try it as a URL
2263
2264    :param str URLorFile: string containing a URL or a file name. Code will try first
2265      to open it as a file and then as a URL.
2266
2267    :returns: a PyCifRW CIF object.
2268    '''
2269    import CifFile as cif # PyCifRW from James Hester
2270
2271    # alternate approach:
2272    #import urllib
2273    #ciffile = 'file:'+urllib.pathname2url(filename)
2274   
2275    try:
2276        fp = open(URLorFile,'r')
2277        cf = cif.ReadCif(fp)
2278        fp.close()
2279        return cf
2280    except IOError:
2281        return cif.ReadCif(URLorFile)
2282
2283if __name__ == '__main__':
2284    app = wx.PySimpleApp()
2285    frm = wx.Frame(None) # create a frame
2286    frm.Show(True)
2287    filename = '/tmp/notzip.zip'
2288    filename = '/tmp/all.zip'
2289    #filename = '/tmp/11bmb_7652.zip'
2290   
2291    #selection=None, confirmoverwrite=True, parent=None
2292    #print ExtractFileFromZip(filename, selection='11bmb_7652.fxye',parent=frm)
2293    print ExtractFileFromZip(filename,multipleselect=True)
2294                             #confirmread=False, confirmoverwrite=False)
2295
2296    # choicelist=[ ('a','b','c'),
2297    #              ('test1','test2'),('no choice',)]
2298    # titles = [ 'a, b or c', 'tests', 'No option here']
2299    # dlg = MultipleChoicesDialog(
2300    #     choicelist,titles,
2301    #     parent=frm)
2302    # if dlg.ShowModal() == wx.ID_OK:
2303    #     print 'Got OK'
Note: See TracBrowser for help on using the repository browser.