source: trunk/imports/G2phase.py @ 3245

Last change on this file since 3245 was 3245, checked in by vondreele, 4 years ago

fix problem with incommensurate magnetic structures from CIFhklReader
make "Edit Body" as "Edit Residue Body" in RB menu
change protein validation bar plots to be vertical
change 'twin' to 'flag' in Reflection List table & add '-2' as 'free' (for proteins - not yet implemented)
fix PDB phase import to get atom type from 76:78 of ATOM record
change importer names for single crystal TOF data to be more explicit (SNS vs ISIS)
change cif structure factor importer - F2, sig(F2) & Fcalc now allowed

  • Property svn:eol-style set to native
  • Property svn:keywords set to Date Author Revision URL Id
File size: 31.3 KB
Line 
1# -*- coding: utf-8 -*-
2########### SVN repository information ###################
3# $Date: 2018-01-27 14:49:12 +0000 (Sat, 27 Jan 2018) $
4# $Author: vondreele $
5# $Revision: 3245 $
6# $URL: trunk/imports/G2phase.py $
7# $Id: G2phase.py 3245 2018-01-27 14:49:12Z vondreele $
8########### SVN repository information ###################
9#
10'''
11*Module G2phase: PDB, .EXP & JANA m40,m50*
12-------------------------------------------
13
14A set of short routines to read in phases using routines that were
15previously implemented in GSAS-II: PDB, GSAS .EXP and JANA m40-m50 file formats
16
17'''
18
19from __future__ import division, print_function
20import sys
21import os.path
22import math
23import random as ran
24import numpy as np
25try:
26    import wx
27except ImportError:
28    wx = None
29import GSASIIobj as G2obj
30import GSASIIspc as G2spc
31import GSASIIlattice as G2lat
32import GSASIIpath
33GSASIIpath.SetVersionNumber("$Revision: 3245 $")
34R2pisq = 1./(2.*np.pi**2)
35
36class PDB_ReaderClass(G2obj.ImportPhase):
37    'Routine to import Phase information from a PDB file'
38    def __init__(self):
39        super(self.__class__,self).__init__( # fancy way to say ImportPhase.__init__
40            extensionlist=('.pdb','.ent','.PDB','.ENT'),
41            strictExtension=True,
42            formatName = 'PDB',
43            longFormatName = 'Original Protein Data Bank (.pdb file) import'
44            )
45    def ContentsValidator(self, filename):
46        '''Taking a stab a validating a PDB file
47        (look for cell & at least one atom)
48        '''
49        fp = open(filename,'r')
50        for i,l in enumerate(fp):
51            if l.startswith('CRYST1'):
52                break
53        else:
54            self.errors = 'no CRYST1 record found'
55            fp.close()
56            return False
57        for i,l in enumerate(fp):
58            if l.startswith('ATOM'):
59                fp.close()
60                return True
61        self.errors = 'no ATOM records found after CRYST1 record'
62        fp.close()
63        return False
64
65    def Reader(self,filename, ParentFrame=None, **unused):
66        'Read a PDF file using :meth:`ReadPDBPhase`'
67        self.Phase = self.ReadPDBPhase(filename, ParentFrame)
68        return True
69       
70    def ReadPDBPhase(self,filename,parent=None):
71        '''Read a phase from a PDB file.
72        '''
73        EightPiSq = 8.*math.pi**2
74        self.errors = 'Error opening file'
75        file = open(filename, 'Ur')
76        Phase = {}
77        Title = ''
78        Compnd = ''
79        Atoms = []
80        A = np.zeros(shape=(3,3))
81        S = file.readline()
82        line = 1
83        SGData = None
84        cell = None
85        while S:
86            self.errors = 'Error reading at line '+str(line)
87            Atom = []
88            if 'TITLE' in S[:5]:
89                Title = S[10:72].strip()
90            elif 'COMPND    ' in S[:10]:
91                Compnd = S[10:72].strip()
92            elif 'CRYST' in S[:5]:
93                abc = S[7:34].split()
94                angles = S[34:55].split()
95                cell=[float(abc[0]),float(abc[1]),float(abc[2]),
96                    float(angles[0]),float(angles[1]),float(angles[2])]
97                Volume = G2lat.calc_V(G2lat.cell2A(cell))
98                AA,AB = G2lat.cell2AB(cell)
99                SpGrp = S[55:65]
100                E,SGData = G2spc.SpcGroup(SpGrp)
101                # space group processing failed, try to look up name in table
102                if E:
103                    SpGrpNorm = G2spc.StandardizeSpcName(SpGrp)
104                    if SpGrpNorm:
105                        E,SGData = G2spc.SpcGroup(SpGrpNorm)
106                while E:
107                    dlg = wx.TextEntryDialog(parent,
108                        SpGrp[:-1]+' is invalid \nN.B.: make sure spaces separate axial fields in symbol',
109                        'ERROR in space group symbol','',style=wx.OK)
110                    if dlg.ShowModal() == wx.ID_OK:
111                        SpGrp = dlg.GetValue()
112                        E,SGData = G2spc.SpcGroup(SpGrp)
113                    else:
114                        SGData = G2obj.P1SGData # P 1
115                        self.warnings += '\nThe space group was not interpreted and has been set to "P 1".'
116                        self.warnings += "Change this in phase's General tab."           
117                    dlg.Destroy()
118#                SGlines = G2spc.SGPrint(SGData)
119#                for l in SGlines: print (l)
120            elif 'SCALE' in S[:5]:
121                V = S[10:41].split()
122                A[int(S[5])-1] = [float(V[0]),float(V[1]),float(V[2])]
123            elif 'ATOM' in S[:4] or 'HETATM' in S[:6]:
124                if not SGData:
125                    self.warnings += '\nThe space group was not read before atoms and has been set to "P 1". '
126                    self.warnings += "Change this in phase's General tab."
127                    SGData = G2obj.P1SGData # P 1
128                XYZ = [float(S[31:39]),float(S[39:47]),float(S[47:55])]
129                XYZ = np.inner(AB,XYZ)
130                XYZ = np.where(abs(XYZ)<0.00001,0,XYZ)
131                SytSym,Mult = G2spc.SytSym(XYZ,SGData)[:2]
132                Uiso = float(S[61:67])/EightPiSq
133                Type = S[76:78].lower()
134                Atom = [S[22:27].strip(),S[17:20].upper(),S[20:22],
135                    S[12:17].strip(),Type.strip().capitalize(),'',XYZ[0],XYZ[1],XYZ[2],
136                    float(S[55:61]),SytSym,Mult,'I',Uiso,0,0,0,0,0,0]
137                if S[16] in [' ','A','B']:      #remove disorered residues - can't handle them just now
138#                    Atom[3] = Atom[3][:3]
139                    Atom.append(ran.randint(0,sys.maxsize))
140                    Atoms.append(Atom)
141            elif 'ANISOU' in S[:6]:
142                Uij = S[30:72].split()
143                Uij = [float(Uij[0])/10000.,float(Uij[1])/10000.,float(Uij[2])/10000.,
144                    float(Uij[3])/10000.,float(Uij[4])/10000.,float(Uij[5])/10000.]
145                Atoms[-1] = Atoms[-1][:14]+Uij
146                Atoms[-1][12] = 'A'
147                Atoms[-1].append(ran.randint(0,sys.maxsize))
148            S = file.readline()
149            line += 1
150        file.close()
151        self.errors = 'Error after read complete'
152        if Title:
153            PhaseName = Title
154        elif Compnd:
155            PhaseName = Compnd
156        else:
157            PhaseName = 'None'
158        if not SGData:
159            raise self.ImportException("No space group (CRYST entry) found")
160        if not cell:
161            raise self.ImportException("No cell (CRYST entry) found")
162        Phase = G2obj.SetNewPhase(Name=PhaseName,SGData=SGData,cell=cell+[Volume,])
163        Phase['General']['Type'] = 'macromolecular'
164        Phase['General']['AtomPtrs'] = [6,4,10,12]
165        Phase['Atoms'] = Atoms
166        return Phase
167
168class EXP_ReaderClass(G2obj.ImportPhase):
169    'Routine to import Phase information from GSAS .EXP files'
170    def __init__(self):
171        super(self.__class__,self).__init__( # fancy way to say ImportPhase.__init__
172            extensionlist=('.EXP','.exp'),
173            strictExtension=True,
174            formatName = 'GSAS .EXP',
175            longFormatName = 'GSAS Experiment (.EXP file) import'
176            )
177       
178    def ContentsValidator(self, filename):
179        'Look for a VERSION tag in 1st line' 
180        fp = open(filename,'r')
181        if fp.read(13) == '     VERSION ':
182            fp.close()
183            return True
184        self.errors = 'File does not begin with VERSION tag'
185        fp.close()
186        return False
187
188    def Reader(self,filename, ParentFrame=None, **unused):
189        'Read a phase from a GSAS .EXP file using :meth:`ReadEXPPhase`'
190        fp = open(filename,'r')
191        self.Phase = self.ReadEXPPhase(ParentFrame, fp)
192        fp.close()
193        return True
194
195    def ReadEXPPhase(self, G2frame,filepointer):
196        '''Read a phase from a GSAS .EXP file.
197        '''
198        shModels = ['cylindrical','none','shear - 2/m','rolling - mmm']
199        textureData = {'Order':0,'Model':'cylindrical','Sample omega':[False,0.0],
200            'Sample chi':[False,0.0],'Sample phi':[False,0.0],'SH Coeff':[False,{}],
201            'SHShow':False,'PFhkl':[0,0,1],'PFxyz':[0,0,1],'PlotType':'Pole figure'}
202        shNcof = 0
203        S = 1
204        NPhas = []
205        Expr = [{},{},{},{},{},{},{},{},{}] # GSAS can have at most 9 phases
206        for line,S in enumerate(filepointer):
207            self.errors = 'Error reading at line '+str(line+1)
208            if 'EXPR NPHAS' in S[:12]:
209                NPhas = S[12:-1].split()
210            if 'CRS' in S[:3]:
211                N = int(S[3:4])-1
212                Expr[N][S[:12]] = S[12:-1]
213        PNames = []
214        if not NPhas:
215            raise self.ImportException("No EXPR NPHAS record found")
216        self.errors = 'Error interpreting file'
217        for n,N in enumerate(NPhas):
218            if N != '0':
219                result = n
220                key = 'CRS'+str(n+1)+'    PNAM'
221                PNames.append(Expr[n][key])
222        if len(PNames) == 0:
223            raise self.ImportException("No phases found")           
224        elif len(PNames) > 1:
225            dlg = wx.SingleChoiceDialog(G2frame, 'Which phase to read?', 'Read phase data', PNames, wx.CHOICEDLG_STYLE)
226            try:
227                if dlg.ShowModal() == wx.ID_OK:
228                    result = dlg.GetSelection() # I think this breaks is there are skipped phases. Cant this happen?
229            finally:
230                dlg.Destroy()       
231        EXPphase = Expr[result]
232        keyList = list(EXPphase.keys())
233        keyList.sort()
234        SGData = {}
235        if NPhas[result] == '1':
236            Ptype = 'nuclear'
237        elif NPhas[result] in ['2','3']:
238            Ptype = 'magnetic'
239            MagDmin = 1.0
240        elif NPhas[result] == '4':
241            Ptype = 'macromolecular'
242        elif NPhas[result] == '10':
243            Ptype = 'Pawley'
244        else:
245            raise self.ImportException("Phase type not recognized")           
246        for key in keyList:
247            if 'PNAM' in key:
248               PhaseName = EXPphase[key].strip()
249            elif 'ABC   ' in key:
250                abc = [float(EXPphase[key][:10]),float(EXPphase[key][10:20]),float(EXPphase[key][20:30])]                       
251            elif 'ANGLES' in key:
252                angles = [float(EXPphase[key][:10]),float(EXPphase[key][10:20]),float(EXPphase[key][20:30])]                                               
253            elif 'SG SYM' in key:
254                SpGrp = EXPphase[key][:15].strip()
255                E,SGData = G2spc.SpcGroup(SpGrp)
256                if E:
257                    SGData = G2obj.P1SGData # P 1 -- unlikely to need this!
258                    self.warnings += '\nThe GSAS space group was not interpreted(!) and has been set to "P 1".'
259                    self.warnings += "Change this in phase's General tab."                       
260            elif 'SPNFLP' in key:
261                SpnFlp = np.array([int(float(s)) for s in EXPphase[key].split()])
262                SpnFlp = np.where(SpnFlp==0,1,SpnFlp)
263                if SGData['SpGrp'][0] in ['A','B','C','I','R','F']:
264                    SpnFlp = list(SpnFlp)+[1,1,1,1]
265            elif 'MXDSTR' in key:
266                MagDmin = float(EXPphase[key][:10])               
267            elif 'OD    ' in key:
268                SHdata = EXPphase[key].split() # may not have all 9 values
269                SHvals = 9*[0]
270                for i in range(9):
271                    try:
272                        float(SHdata[i])
273                        SHvals[i] = SHdata[i]
274                    except:
275                        pass
276                textureData['Order'] = int(SHvals[0])
277                textureData['Model'] = shModels[int(SHvals[2])]
278                textureData['Sample omega'] = [False,float(SHvals[6])]
279                textureData['Sample chi'] = [False,float(SHvals[7])]
280                textureData['Sample phi'] = [False,float(SHvals[8])]
281                shNcof = int(SHvals[1])
282        Atoms = []
283        Amat,Bmat = G2lat.cell2AB(abc+angles)
284        if Ptype in ['nuclear','magnetic',]:
285            for key in keyList:
286                if 'AT' in key:
287                    if key[11:] == 'A':
288                        S = EXPphase[key]
289                    elif key[11:] == 'B':
290                        S1 = EXPphase[key]
291                        Atom = [S[50:58].strip(),S[:10].strip().capitalize(),'',
292                            float(S[10:20]),float(S[20:30]),float(S[30:40]),
293                            float(S[40:50]),'',int(S[60:62]),S1[62:63]]
294                            #float(S[40:50]),'',int(S[60:62]),S1[130:131]]
295                        if Atom[9] == 'I':
296                            Atom += [float(S1[0:10]),0.,0.,0.,0.,0.,0.]
297                        elif Atom[9] == 'A':
298                            Atom += [0.0,
299                                float(S1[ 0:10]),float(S1[10:20]),
300                                float(S1[20:30]),float(S1[30:40]),
301                                float(S1[40:50]),float(S1[50:60])]
302                        else:
303                            print('Error in line with key: '+key)
304                            Atom += [0.,0.,0.,0.,0.,0.,0.]
305                        XYZ = Atom[3:6]
306                        Atom[7],Atom[8] = G2spc.SytSym(XYZ,SGData)[:2]
307                        Atom.append(ran.randint(0,sys.maxsize))
308                        Atoms.append(Atom)
309                    elif key[11:] == 'M' and key[6:8] == 'AT':
310                        Ptype = 'magnetic'
311                        S = EXPphase[key]
312                        mom = np.array([float(S[:10]),float(S[10:20]),float(S[20:30])])
313                        mag = np.sqrt(np.sum(mom**2))
314                        mom = np.inner(Bmat,mom)*mag
315                        Atoms[-1] = Atom[:7]+list(mom)+Atom[7:]
316        elif Ptype == 'macromolecular':
317            for key in keyList:
318                if 'AT' in key[6:8]:
319                    S = EXPphase[key]
320                    Atom = [S[56:60].strip(),S[50:54].strip().upper(),S[54:56],
321                        S[46:51].strip(),S[:8].strip().capitalize(),'',
322                        float(S[16:24]),float(S[24:32]),float(S[32:40]),
323                        float(S[8:16]),'1',1,'I',float(S[40:46]),0,0,0,0,0,0]
324                    XYZ = Atom[6:9]
325                    Atom[10],Atom[11] = G2spc.SytSym(XYZ,SGData)[:2]
326                    Atom.append(ran.randint(0,sys.maxsize))
327                    Atoms.append(Atom)
328        Volume = G2lat.calc_V(G2lat.cell2A(abc+angles))
329        if shNcof:
330            shCoef = {}
331            nRec = [i+1 for i in range((shNcof-1)/6+1)]
332            for irec in nRec:
333                ODkey = keyList[0][:6]+'OD'+'%3dA'%(irec)
334                indx = EXPphase[ODkey].split()
335                ODkey = ODkey[:-1]+'B'
336                vals = EXPphase[ODkey].split()
337                for i,val in enumerate(vals):
338                    key = 'C(%s,%s,%s)'%(indx[3*i],indx[3*i+1],indx[3*i+2])
339                    shCoef[key] = float(val)
340            textureData['SH Coeff'] = [False,shCoef]
341        if not SGData:
342            raise self.ImportException("No space group found in phase")
343        if not abc:
344            raise self.ImportException("No cell lengths found in phase")
345        if not angles:
346            raise self.ImportException("No cell angles found in phase")
347        if not Atoms:
348            raise self.ImportException("No atoms found in phase")
349        Phase = G2obj.SetNewPhase(Name=PhaseName,SGData=SGData,cell=abc+angles+[Volume,])
350        general = Phase['General']
351        general['Type'] = Ptype
352        if general['Type'] =='macromolecular':
353            general['AtomPtrs'] = [6,4,10,12]
354        elif general['Type'] =='magnetic':
355            general['AtomPtrs'] = [3,1,10,12]
356            general['SGData']['SGSpin'] = SpnFlp
357            general['MagDmin'] = MagDmin   
358        else:   #nuclear
359            general['AtomPtrs'] = [3,1,7,9]   
360        general['SH Texture'] = textureData
361        Phase['Atoms'] = Atoms
362        return Phase
363
364class JANA_ReaderClass(G2obj.ImportPhase):
365    'Routine to import Phase information from a JANA2006 file'
366    def __init__(self):
367        super(self.__class__,self).__init__( # fancy way to say ImportPhase.__init__
368            extensionlist=('.m50','.M50'),
369            strictExtension=True,
370            formatName = 'JANA m50',
371            longFormatName = 'JANA2006 phase import'
372            )
373    def ContentsValidator(self, filename):
374        '''Taking a stab a validating a .m50 file
375        (look for cell & at least one atom)
376        '''
377        fp = open(filename,'r')
378        for i,l in enumerate(fp):
379            if l.startswith('cell'):
380                break
381        else:
382            self.errors = 'no cell record found'
383            fp.close()
384            return False
385        for i,l in enumerate(fp):
386            if l.startswith('spgroup'):
387                fp.close()
388                return True
389        self.errors = 'no spgroup record found after cell record'
390        fp.close()
391        return False
392       
393    def Reader(self,filename, ParentFrame=None, **unused):
394        'Read a m50 file using :meth:`ReadJANAPhase`'
395        self.Phase = self.ReadJANAPhase(filename, ParentFrame)
396        return True
397       
398    def ReadJANAPhase(self,filename,parent=None):
399        '''Read a phase from a JANA2006 m50 & m40 files.
400        '''
401        self.errors = 'Error opening file'
402        fp = open(filename, 'Ur') #contains only cell & spcgroup
403        Phase = {}
404        Title = os.path.basename(filename)
405        Type = 'nuclear'
406        Atoms = []
407        Atypes = []
408        SuperVec = [[0,0,.1],False,4]
409        S = fp.readline()
410        line = 1
411        SGData = None
412        SuperSg = ''
413        cell = None
414        nqi = 0
415        version = '2000'
416        while S:
417            self.errors = 'Error reading at line '+str(line)
418            if 'title' in S and S != 'title\n':
419                Title = S.split()[1]
420            elif 'Jana2006' in S:
421                self.warnings += '\nJana2006 file detected'
422                version = '2006'
423            elif 'cell' in S[:4]:
424                cell = S[5:].split()
425                cell=[float(cell[0]),float(cell[1]),float(cell[2]),
426                    float(cell[3]),float(cell[4]),float(cell[5])]
427                Volume = G2lat.calc_V(G2lat.cell2A(cell))
428                G,g = G2lat.cell2Gmat(cell)
429                ast = np.sqrt(np.diag(G))
430                Mast = np.multiply.outer(ast,ast)   
431               
432            elif 'spgroup' in S:
433                if 'X' in S:
434                    raise self.ImportException("Ad hoc Supersymmetry centering "+S+" not allowed in GSAS-II")           
435                SpGrp = S.split()[1]
436                SuperSg = ''
437                if '(' in SpGrp:    #supercell symmetry - split in 2
438                    SuperStr = SpGrp.split('(')
439                    SpGrp = SuperStr[0]
440                    SuperSg = '('+SuperStr[1]
441                SpGrpNorm = G2spc.StandardizeSpcName(SpGrp)
442                E,SGData = G2spc.SpcGroup(SpGrpNorm)
443                # space group processing failed, try to look up name in table
444                while E:
445                    print (G2spc.SGErrors(E))
446                    dlg = wx.TextEntryDialog(parent,
447                        SpGrp[:-1]+' is invalid \nN.B.: make sure spaces separate axial fields in symbol',
448                        'ERROR in space group symbol','',style=wx.OK)
449                    if dlg.ShowModal() == wx.ID_OK:
450                        SpGrp = dlg.GetValue()
451                        E,SGData = G2spc.SpcGroup(SpGrp)
452                    else:
453                        SGData = G2obj.P1SGData # P 1
454                        self.warnings += '\nThe space group was not interpreted and has been set to "P 1".'
455                        self.warnings += "Change this in phase's General tab."           
456                    dlg.Destroy()
457                G2spc.SGPrint(SGData) #silent check of space group symbol
458            elif 'qi' in S[:2]:
459                if nqi:
460                    raise self.ImportException("Supersymmetry too high; GSAS-II limited to (3+1) supersymmetry")           
461                vec = S.split()[1:]
462                SuperVec = [[float(vec[i]) for i in range(3)],False,4]
463                nqi += 1
464            elif 'atom' in S[:4]:
465                Atypes.append(S.split()[1])
466            S = fp.readline()
467            line += 1
468        fp.close()
469        #read atoms from m40 file
470        if not SGData:
471            self.warnings += '\nThe space group was not read before atoms and has been set to "P 1". '
472            self.warnings += "Change this in phase's General tab."
473            SGData = G2obj.P1SGData # P 1
474        waveTypes = ['Fourier','Sawtooth','ZigZag',]
475        filename2 = os.path.splitext(filename)[0]+'.m40'
476        file2 = open(filename2,'Ur')
477        S = file2.readline()
478        line = 1
479        self.errors = 'Error reading at line '+str(line)
480        nAtoms = int(S.split()[0])
481        for i in range(4):
482            S = file2.readline()           
483        for i in range(nAtoms):
484            S1 = file2.readline()
485            S1N = S1.split()[-3:]   # no. occ, no. pos waves, no. ADP waves
486            S1N = [int(i) for i in S1N]
487            S1T = list(S1[60:63])
488            waveType = waveTypes[int(S1T[1])]
489            Spos = []
490            Sadp = []
491            Sfrac = []
492            Smag = []
493            XYZ = [float(S1[27:36]),float(S1[36:45]),float(S1[45:54])]
494            SytSym,Mult = G2spc.SytSym(XYZ,SGData)[:2]
495            aType = Atypes[int(S1[9:11])-1]
496            Name = S1[:8].strip()
497            if S1[11:15].strip() == '1':
498                S2 = file2.readline()
499                Uiso = S2[:9]
500                if version == '2000':
501                    Uiso = R2pisq*float(Uiso)/4.      #Biso -> Uiso
502                Uij = [0,0,0,0,0,0]
503                IA = 'I'
504            elif S1[11:15].strip() == '2':
505                S2 = file2.readline()
506                IA = 'A'
507                Uiso = 0.
508                Uij = [float(S2[:9]),float(S2[9:18]),float(S2[18:27]),
509                    float(S2[27:36]),float(S2[36:45]),float(S2[45:54])] #Uij in Jana2006!
510                if version == '2000':
511                    Uij = R2pisq*G2lat.UijtoU6(G2lat.U6toUij(Uij)/Mast) #these things are betaij in Jana2000! need to convert to Uij
512            for i in range(S1N[0]):
513                if not i:
514                    FS = file2.readline()
515                    Sfrac.append(FS[:9])    #'O' or 'delta' = 'length' for crenel
516                    if int(S1T[0]):  #"", "Legendre" or "Xharm" in 18:27 for "crenel"!
517                        waveType = 'Crenel/Fourier' #all waves 'Fourier' no other choice
518                Sfrac.append(file2.readline()[:18]) #if not crenel = Osin & Ocos
519                # else Osin & Ocos except last one is X40 = 'Center'
520            for i in range(S1N[1]): 
521                Spos.append(file2.readline()[:54])
522            for i in range(S1N[2]):
523                Sadp.append(file2.readline()[:54]+file2.readline())
524            if sum(S1N):    #if any waves: skip mystery line?
525                file2.readline()
526            for i,it in enumerate(Sfrac):
527                print (i,it)
528                if not i:
529                    if 'Crenel' in waveType:
530                        vals = [float(it),float(Sfrac[-1][:9])]
531                    else:
532                        vals = [float(it),]
533                else:
534                    vals = [float(it[:9]),float(it[9:18])]
535                if 'Crenel' in waveType and i == len(Sfrac)-1:
536                    del Sfrac[-1]
537                    break               
538                Sfrac[i] = [vals,False]
539                print (Sfrac[i])
540            for i,it in enumerate(Spos):
541                if waveType in ['Sawtooth',] and not i:
542                    vals = [float(it[:9]),float(it[9:18]),float(it[18:27]),float(it[27:36])]
543                else:
544                    vals = [float(it[:9]),float(it[9:18]),float(it[18:27]),float(it[27:36]),float(it[36:45]),float(it[45:54])]
545                Spos[i] = [vals,False]
546            for i,it in enumerate(Sadp):
547                vals = [float(it[:9]),float(it[9:18]),float(it[18:27]),float(it[27:36]),float(it[36:45]),float(it[45:54]),
548                    float(it[54:63]),float(it[63:72]),float(it[72:81]),float(it[81:90]),float(it[90:99]),float(it[99:108])]
549                #these are betaij modulations in Jana2000! need to convert to Uij modulations
550                if version == '2000':               
551                    vals[:6] = R2pisq*G2lat.UijtoU6(G2lat.U6toUij(vals[:6])/Mast)    #convert sin bij to Uij
552                    vals[6:] = R2pisq*G2lat.UijtoU6(G2lat.U6toUij(vals[6:])/Mast)    #convert cos bij to Uij
553                Sadp[i] = [vals,False]
554            Atom = [Name,aType,'',XYZ[0],XYZ[1],XYZ[2],1.0,SytSym,Mult,IA,Uiso]
555            Atom += Uij
556            Atom.append(ran.randint(0,sys.maxsize))
557            Atom.append([])
558            Atom.append([])
559            Atom.append({'SS1':{'waveType':waveType,'Sfrac':Sfrac,'Spos':Spos,'Sadp':Sadp,'Smag':Smag}})    #SS2 is for (3+2), etc.
560            Atoms.append(Atom)
561        file2.close()
562        self.errors = 'Error after read complete'
563        if not SGData:
564            raise self.ImportException("No space group (spcgroup entry) found")
565        if not cell:
566            raise self.ImportException("No cell found")
567        Phase = G2obj.SetNewPhase(Name=Title,SGData=SGData,cell=cell+[Volume,])
568        Phase['General']['Type'] = Type
569        Phase['General']['Modulated'] = True
570        Phase['General']['Super'] = nqi
571        Phase['General']['SuperVec'] = SuperVec
572        Phase['General']['SuperSg'] = SuperSg
573        if SuperSg:
574            Phase['General']['SSGData'] = G2spc.SSpcGroup(SGData,SuperSg)[1]
575        Phase['General']['AtomPtrs'] = [3,1,7,9]
576        Phase['Atoms'] = Atoms
577        return Phase
578   
579class PDF_ReaderClass(G2obj.ImportPhase):
580    'Routine to import Phase information from ICDD PDF Card files'
581    def __init__(self):
582        super(self.__class__,self).__init__( # fancy way to say ImportPhase.__init__
583            extensionlist=('.str',),
584            strictExtension=True,
585            formatName = 'ICDD .str',
586            longFormatName = 'ICDD PDF Card (.str file) import'
587            )
588       
589    def ContentsValidator(self, filename):
590        'Look for a str tag in 1st line' 
591        fp = open(filename,'r')
592        if fp.read(3) == 'str':
593            fp.close()
594            return True
595        self.errors = 'File does not begin with str tag'
596        fp.close()
597        return False
598
599    def Reader(self,filename, ParentFrame=None, **unused):
600        'Read phase from a ICDD .str file using :meth:`ReadPDFPhase`'
601        fp = open(filename,'r')
602        self.Phase = self.ReadPDFPhase(ParentFrame, fp)
603        fp.close()
604        return True
605
606    def ReadPDFPhase(self, G2frame,fp):
607        '''Read a phase from a ICDD .str file.
608        '''
609        EightPiSq = 8.*math.pi**2
610        self.errors = 'Error opening file'
611        Phase = {}
612        Atoms = []
613        S = fp.readline()
614        line = 1
615        SGData = None
616        cell = []
617        cellkey = []
618        while S:
619            if 'space_group' in S:
620                break
621            S = fp.readline()                   
622        while S:
623            self.errors = 'Error reading at line '+str(line)
624            if 'phase_name' in S:
625                Title = S.split('"')[1]
626            elif 'Space group (HMS)' in S:
627                SpGrp = S.split()[-1]
628                SpGrpNorm = G2spc.StandardizeSpcName(SpGrp)
629                E,SGData = G2spc.SpcGroup(SpGrpNorm)
630                # space group processing failed, try to look up name in table
631                while E:
632                    print (G2spc.SGErrors(E))
633                    dlg = wx.TextEntryDialog(G2frame,
634                        SpGrp[:-1]+' is invalid \nN.B.: make sure spaces separate axial fields in symbol',
635                        'ERROR in space group symbol','',style=wx.OK)
636                    if dlg.ShowModal() == wx.ID_OK:
637                        SpGrp = dlg.GetValue()
638                        E,SGData = G2spc.SpcGroup(SpGrp)
639                    else:
640                        SGData = G2obj.P1SGData # P 1
641                        self.warnings += '\nThe space group was not interpreted and has been set to "P 1".'
642                        self.warnings += "Change this in phase's General tab."           
643                    dlg.Destroy()
644                G2spc.SGPrint(SGData) #silent check of space group symbol
645            elif 'a a_' in S[:7]:
646                data = S.split()
647                cell.append(float(data[2]))
648                cellkey.append(data[1])
649            elif 'b b_' in S[:7]:
650                data = S.split()
651                cell.append(float(data[2]))
652                cellkey.append(data[1])
653            elif 'b =' in S[:6]:
654                data = S.split('=')
655                indx = cellkey.index(data[1].split(';')[0])
656                cell.append(cell[indx])
657            elif 'c c_' in S[:7]:
658                data = S.split()
659                cell.append(float(data[2]))
660            elif 'c =' in S[:6]:
661                data = S.split('=')
662                indx = cellkey.index(data[1].split(';')[0])
663                cell.append(cell[indx])
664            elif 'al' in S[:5]:
665                cell.append(float(S.split()[1]))
666            elif 'be' in S[:5]:
667                cell.append(float(S.split()[1]))
668            elif 'ga' in S[:5]:
669                cell.append(float(S.split()[1]))
670                Volume = G2lat.calc_V(G2lat.cell2A(cell))
671                break
672            S = fp.readline()
673        S = fp.readline()
674        while S:
675            if '/*' in S[:5]:
676                break
677            if 'site' in S[:7]:
678                atom = []
679                xyzkey = []
680                data = S.split()
681                atom.append(data[1])    #name
682                pos = data.index('occ')+1
683                atom.append(data[pos])  #type
684                atom.append('')         #refine
685                for xid in ['x =','y =','z =']:
686                    if xid in S:
687                        xpos = S.index(xid)+3
688                        xend = xpos+S[xpos:].index(';')
689                        if S[xpos:xend] in xyzkey:
690                            atom.append(atom[3+xyzkey.index(S[xpos:xend])])
691                        else:
692                            atom.append(eval(S[xpos:xend]+'.'))
693                    else:
694                        xpos = data.index(xid[0])+2
695                        xyzkey.append(data[xpos-1][1:])
696                        atom.append(float(data[xpos]))
697                atom.append(float(data[pos+2]))
698                SytSym,Mult = G2spc.SytSym(np.array(atom[3:6]),SGData)[:2]
699                atom.append(SytSym)
700                atom.append(Mult)
701                if 'beq' in S:
702                    atom.append('I')
703                    upos = data.index('beq')
704                    atom.append(float(data[upos+2])/EightPiSq)
705                    atom += [0.,0.,0.,0.,0.,0.,]
706                elif 'ADPs' in S:
707                    upos = data.index('ADPs')
708                    atom.append('A')
709                    atom.append(0.0)
710                    for uid in ['Bani11','Bani22','Bani33','Bani12','Bani13','Bani23']:
711                        upos = data.index(uid)+1
712                        atom.append(float(data[upos])/EightPiSq)
713                else:
714                    atom.append('I')
715                    atom += [0.02,0.,0.,0.,0.,0.,0.,]                   
716                atom.append(ran.randint(0,sys.maxsize))
717                Atoms.append(atom)
718            S = fp.readline()               
719        fp.close()
720        self.errors = 'Error after read complete'
721        if not SGData:
722            raise self.ImportException("No space group (spcgroup entry) found")
723        if not cell:
724            raise self.ImportException("No cell found")
725        Phase = G2obj.SetNewPhase(Name=Title,SGData=SGData,cell=cell+[Volume,])
726        Phase['General']['Type'] = 'nuclear'
727        Phase['General']['AtomPtrs'] = [3,1,7,9]
728        Phase['Atoms'] = Atoms
729        return Phase
Note: See TracBrowser for help on using the repository browser.