source: trunk/imports/G2phase.py @ 3221

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

fix import of magnetic structures from .EXP files
complete list of pre & post 2002 cubic space groups
make sure SGDataSGGray? is present

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