source: trunk/GSASIIspc.py @ 342

Last change on this file since 342 was 342, checked in by vondreele, 12 years ago

major modifications to get 1st "working" version of Refine
Add GSASIImapvars.py
fix cell refinement in indexing window
single cycle for peak refinement
  • Property svn:keywords set to Date Author Revision URL Id
File size: 54.9 KB
Line 
1"GSASII - Space group interpretion routines"
2
3########### SVN repository information ###################
4# $Date: 2011-08-05 19:35:43 +0000 (Fri, 05 Aug 2011) $
5# $Author: vondreele $
6# $Revision: 342 $
7# $URL: trunk/GSASIIspc.py $
8# $Id: GSASIIspc.py 342 2011-08-05 19:35:43Z vondreele $
9########### SVN repository information ###################
10import numpy as np
11import numpy.ma as ma
12import numpy.linalg as nl
13import math
14import sys
15import os.path as ospath
16
17import GSASIIpath
18import pyspg
19
20def SpcGroup(SGSymbol):
21    '''
22    Determines cell and symmetry information from a short H-M space group name
23    input:
24        SGSymbol - space group symbol (string) with spaces between axial fields
25    returns:
26        SGError = 0 for no errors; >0 for errors (see SGErrors below for details)
27        SGData - dictionary with entries:
28             'SpGrp': space group symbol slightly cleaned up
29             'Laue':  one of '-1','2/m','mmm','4/m','4/mmm','3R','3mR','3',
30                      '3m1','31m','6/m','6/mmm','m3','m3m'
31             'SGInv': boolean; True if centrosymmetric, False if not
32             'SGLatt': one of 'P','A','B','C','I','F','R'
33             'SGUniq': one of 'a','b','c' if monoclinic, '' otherwise
34             'SGCen': cell centering vectors [0,0,0] at least
35             'SGOps': symmetry operations as [M,T] so that M*x+T = x'
36             'SGSys': one of 'triclinic','monoclinic','orthorhombic','tetragonal','rhombohedral','trigonal','hexagonal','cubic'
37       '''
38    LaueSym = ('-1','2/m','mmm','4/m','4/mmm','3R','3mR','3','3m1','31m','6/m','6/mmm','m3','m3m')
39    LattSym = ('P','A','B','C','I','F','R')
40    UniqSym = ('','','a','b','c','',)
41    SysSym = ('triclinic','monoclinic','orthorhombic','tetragonal','rhombohedral','trigonal','hexagonal','cubic')
42    SGData = {}
43    SGData['SpGrp'] = SGSymbol.strip().lower().capitalize()
44    SGInfo = pyspg.sgforpy(SGSymbol)
45    SGData['SGLaue'] = LaueSym[SGInfo[0]-1]
46    SGData['SGInv'] = bool(SGInfo[1])
47    SGData['SGLatt'] = LattSym[SGInfo[2]-1]
48    SGData['SGUniq'] = UniqSym[SGInfo[3]+1]
49    if SGData['SGLatt'] == 'P':
50        SGData['SGCen'] = np.array(([0,0,0],))
51    elif SGData['SGLatt'] == 'A':
52        SGData['SGCen'] = np.array(([0,0,0],[0,.5,.5]))
53    elif SGData['SGLatt'] == 'B':
54        SGData['SGCen'] = np.array(([0,0,0],[.5,0,.5]))
55    elif SGData['SGLatt'] == 'C':
56        SGData['SGCen'] = np.array(([0,0,0],[.5,.5,0,]))
57    elif SGData['SGLatt'] == 'I':
58        SGData['SGCen'] = np.array(([0,0,0],[.5,.5,.5]))
59    elif SGData['SGLatt'] == 'F':
60        SGData['SGCen'] = np.array(([0,0,0],[0,.5,.5],[.5,0,.5],[.5,.5,0,]))
61    elif SGData['SGLatt'] == 'R':
62        SGData['SGCen'] = np.array(([0,0,0],[1./3.,2./3.,2./3.],[2./3.,1./3.,1./3.]))
63    SGData['SGOps'] = []
64    for i in range(SGInfo[5]):
65        Mat = np.array(SGInfo[6][i])
66        Trns = np.array(SGInfo[7][i])
67        SGData['SGOps'].append([Mat,Trns])
68    if SGData['SGLaue'] in '-1':
69        SGData['SGSys'] = SysSym[0]
70    elif SGData['SGLaue'] in '2/m':
71        SGData['SGSys'] = SysSym[1]
72    elif SGData['SGLaue'] in 'mmm':
73        SGData['SGSys'] = SysSym[2]
74    elif SGData['SGLaue'] in ['4/m','4/mmm']:
75        SGData['SGSys'] = SysSym[3]
76    elif SGData['SGLaue'] in ['3R','3mR']:
77        SGData['SGSys'] = SysSym[4]
78    elif SGData['SGLaue'] in ['3','3m1','31m']:
79        SGData['SGSys'] = SysSym[5]
80    elif SGData['SGLaue'] in ['6/m','6/mmm']:
81        SGData['SGSys'] = SysSym[6]
82    elif SGData['SGLaue'] in ['m3','m3m']:
83        SGData['SGSys'] = SysSym[7]
84    return SGInfo[8],SGData
85
86def SGErrors(IErr):
87    '''
88    Interprets the error message code from SpcGroup. Used in SpaceGroup.
89    input:
90        SGError - from SpcGroup
91    returns:
92        ErrString - a string with the error message or "Unknown error"
93    '''
94
95    ErrString = [' ',
96        'Less than 2 operator fields were found',
97        'Illegal Lattice type, not P, A, B, C, I, F or R',
98        'Rhombohedral lattice requires a 3-axis',
99        'Minus sign does not preceed 1, 2, 3, 4 or 6',
100        'Either a 5-axis anywhere or a 3-axis in field not allowed',
101        ' ',
102        'I for COMPUTED GO TO out of range.',
103        'An a-glide mirror normal to A not allowed',
104        'A b-glide mirror normal to B not allowed',
105        'A c-glide mirror normal to C not allowed',
106        'D-glide in a primitive lattice not allowed',
107        'A 4-axis not allowed in the 2nd operator field',
108        'A 6-axis not allowed in the 2nd operator field',
109        'More than 24 matrices needed to define group',
110        ' ',
111        'Improper construction of a rotation operator',
112        'Mirror following a / not allowed',
113        'A translation conflict between operators',
114        'The 2bar operator is not allowed',
115        '3 fields are legal only in R & m3 cubic groups',
116        'Syntax error. Expected I -4 3 d at this point',
117        ' ',
118        'A or B centered tetragonal not allowed',
119        ' ','unknown error in sgroup',' ',' ',' ',
120        'Illegal character in the space group symbol',
121        ]
122    try:
123        return ErrString[IErr]
124    except:
125        return "Unknown error"
126   
127def SGPrint(SGData):
128    '''
129    Print the output of SpcGroup in a nicely formatted way. Used in SpaceGroup
130    input:
131        SGData - from SpcGroup
132    returns:
133        SGText - list of strings with the space group details
134    '''
135    POL = (' ','x','y','x y','z','x z','y z','xyz','111')
136    Mult = len(SGData['SGCen'])*len(SGData['SGOps'])*(int(SGData['SGInv'])+1)
137    NP = [1,2,4]
138    NPZ = [0,1]
139    for M,T in SGData['SGOps']:
140        for i in range(3):
141            if M[i][i] <= 0.: NP[i] = 0
142        if M[0][2] > 0: NPZ[0] = 8
143        if M[1][2] > 0: NPZ[1] = 0
144    NPol = (NP[0]+NP[1]+NP[2]+NPZ[0]*NPZ[1])*(1-int(SGData['SGInv']))
145    SGText = []
146    SGText.append(' Space Group: '+SGData['SpGrp'])
147    CentStr = 'centrosymmetric'
148    if not SGData['SGInv']:
149        CentStr = 'non'+CentStr
150    if SGData['SGLatt'] in 'ABCIFR':
151        SGText.append(' The lattice is '+CentStr+' '+SGData['SGLatt']+'-centered '+SGData['SGSys'].lower())
152    else:
153        SGText.append(' The lattice is '+CentStr+' '+'primitive '+SGData['SGSys'].lower())       
154    SGText.append(' Multiplicity of a general site is '+str(Mult))
155    SGText.append(' The Laue symmetry is '+SGData['SGLaue'])
156    if SGData['SGUniq'] in ['a','b','c']:
157        SGText.append(' The unique monoclinic axis is '+SGData['SGUniq'])
158    if SGData['SGInv']:
159        SGText.append(' The inversion center is located at 0,0,0')
160    if NPol:
161        SGText.append(' The location of the origin is arbitrary in '+POL[NPol])
162    SGText.append('\n'+' The equivalent positions are:')
163    if SGData['SGLatt'] != 'P':
164        SGText.append('\n ('+Latt2text(SGData['SGLatt'])+')+')
165    if SGData['SGLaue'] in ['-1','2/m','mmm','4/m','4/mmm']:
166        Ncol = 2
167    else:
168        Ncol = 3
169    line = ' '
170    for iop,[M,T] in enumerate(SGData['SGOps']):
171        if iop % Ncol == 0:
172            SGText.append(line)       
173            line = ' '
174        Fld = '(%2i) ' % (iop+1)+MT2text(M,T)+'\t'
175        line += Fld
176    SGText.append(line)       
177    return SGText
178   
179def MT2text(M,T):
180    #From space group matrix/translation operator returns text version
181    XYZ = ('-Z ','-Y ','-X ','X-Y','ERR','Y-X',' X ',' Y ',' Z ','+X ','+Y ','+Z ')
182    TRA = ('   ','ERR','1/6','1/4','1/3','ERR','1/2','ERR','2/3','3/4','5/6','ERR')
183    Fld = ''
184    for j in range(3):
185        IJ = int(round(2*M[j][0]+3*M[j][1]+4*M[j][2]+4)) % 12
186        IK = int(round(T[j]*12)) % 12
187        if IK > 0 and IJ > 4: IJ += 3
188        Fld += TRA[IK]+XYZ[IJ]
189        if j != 2: Fld += ','
190    return Fld
191   
192def Latt2text(Latt):
193    #From lattice type ('P',A', etc.) returns ';' delimited cell centering vectors
194    lattTxt = {'A':'0,0,0; 0,1/2,1/2','B':'0,0,0; 1/2,0,1/2',
195        'C':'0,0,0; 1/2,1/2,0','I':'0,0,0; 1/2,1/2,1/2',
196        'F':'0,0,0; 0,1/2,1/2; 1/2,0,1/2; 1/2,1/2,0',
197        'R':'0,0,0; 1/3,2/3,2/3; 2/3,1/3,1/3','P':'0,0,0'}
198    return lattTxt[Latt]   
199       
200def SpaceGroup(SGSymbol):
201    '''
202    Print the output of SpcGroup in a nicely formatted way.
203    input:
204        SGSymbol - space group symbol (string) with spaces between axial fields
205    returns:
206        nothing
207    '''
208    E,A = SpcGroup(SGSymbol)
209    if E > 0:
210        print SGErrors(E)
211        return
212    for l in SGPrint(A):
213        print l
214
215def MoveToUnitCell(xyz):
216    '''
217    Translates a set of coordinates so that all values are >=0 and < 1
218    input:
219        xyz - a list or numpy array of fractional coordinates
220    returns:
221        XYZ - numpy array of new coordinates inside 0-1
222    '''
223    XYZ = np.zeros(3)
224    for i,x in enumerate(xyz):
225        XYZ[i] = (x-int(x))%1.0
226    return XYZ
227       
228def Opposite(XYZ,toler=0.0002):
229    '''
230    Gives opposite corner, edge or face of unit cell for position within tolerance.
231        Result may be just outside the cell within tolerance
232    input:
233        XYZ: 0 >= np.array[x,y,z] > 1 as by MoveToUnitCell
234        toler: unit cell fraction tolerance making opposite
235    returns:
236        XYZ: array of opposite positions; always contains XYZ
237    '''
238    perm3 = [[1,1,1],[0,1,1],[1,0,1],[1,1,0],[1,0,0],[0,1,0],[0,0,1],[0,0,0]]
239    TB = np.where(abs(XYZ-1)<toler,-1,0)+np.where(abs(XYZ)<toler,1,0)
240    perm = TB*perm3
241    cperm = ['%d%d%d'%(i,j,k) for i,j,k in perm]
242    D = dict(zip(cperm,perm))
243    new = []
244    for key in D:
245        new.append(np.array(D[key])+np.array(XYZ))
246    return new
247       
248def GenAtom(XYZ,SGData,All=False,Uij=[]):
249    '''
250    Generates the equivalent positions for a specified coordinate and space group
251    input: 
252        XYZ an array, tuple or list containing 3 elements: x, y & z
253        SGData, from SpcGroup
254        All  = True return all equivalent positions including duplicates
255             = False return only unique positions
256        Uij  = [U11,U22,U33,U12,U13,U23] or [] if no Uij       
257    return: [[XYZEquiv],Idup,[UijEquiv]]
258        [XYZEquiv] is list of equivalent positions (XYZ is first entry)
259        Idup = [-][C]SS where SS is the symmetry operator number (1-24), C (if not 0,0,0)
260        is centering operator number (1-4) and - is for inversion
261        Cell = unit cell translations needed to put new positions inside cell
262        [UijEquiv] - equivalent Uij; absent if no Uij given
263    '''
264    XYZEquiv = []
265    UijEquiv = []
266    Idup = []
267    Cell = []
268    X = np.array(XYZ)
269    X = MoveToUnitCell(X)
270    for ic,cen in enumerate(SGData['SGCen']):
271        C = np.array(cen)
272        for invers in range(int(SGData['SGInv']+1)):
273            for io,[M,T] in enumerate(SGData['SGOps']):
274                idup = ((io+1)+100*ic)*(1-2*invers)
275                XT = np.inner(M,X)+T
276                if len(Uij):
277                    U = Uij2U(Uij)
278                    U = np.inner(M,np.inner(U,M).T)
279                    newUij = U2Uij(U)
280                if invers:
281                    XT = -XT
282                XT += C
283                newX = MoveToUnitCell(XT)
284                cell = np.asarray(np.rint(newX-XT),dtype=np.int32)
285                if All:
286                    if np.allclose(newX,X,atol=0.0002):
287                        idup = False
288                else:
289                    if True in [np.allclose(newX,X,atol=0.0002) for oldX in XYZEquiv]:
290                        idup = False
291                if All or idup:
292                    XYZEquiv.append(newX)
293                    Idup.append(idup)
294                    Cell.append(cell)
295                    if len(Uij):
296                        UijEquiv.append(newUij)                   
297    if len(Uij):
298        return zip(XYZEquiv,UijEquiv,Idup,Cell)
299    else:
300        return zip(XYZEquiv,Idup,Cell)
301
302def GenHKLf(HKL,SGData,Friedel=False):
303    '''
304    Uses old GSAS Fortran routine genhkl.for
305    input:
306        HKL - [h,k,l]
307        SGData - space group data obtained from SpcGroup
308        Friedel = True to retain Friedel pairs for centrosymmetric case
309    returns:
310        iabsnt = True is reflection is forbidden by symmetry
311        mulp = reflection multiplicity including Fridel pairs
312        Uniq = numpy array of equivalent hkl in descending order of h,k,l
313    '''
314    hklf = HKL+[0,]
315    Ops = SGData['SGOps']
316    OpM = np.array([op[0] for op in Ops])
317    OpT = np.array([op[1] for op in Ops])
318    Inv = SGData['SGInv']
319    Cen = np.array([cen for cen in SGData['SGCen']])
320   
321    Nuniq,Uniq,iabsnt,mulp = pyspg.genhklpy(hklf,len(Ops),OpM,OpT,SGData['SGInv'],len(Cen),Cen)
322    h,k,l,f = Uniq
323    Uniq=np.array(zip(h[:Nuniq],k[:Nuniq],l[:Nuniq]))
324    Uniq = np.array(Uniq)
325   
326    return iabsnt,2*mulp,Uniq       #include Friedel pairs in powder mulp
327
328       
329def GenHKL(HKL,SGData,Friedel=False):          #gives wrong answers!
330    '''
331    Generate set of equivalent reflections
332    input:
333        HKL - [h,k,l]
334        SGData - space group data obtained from SpcGroup
335        Friedel = True to retain Friedel pairs for centrosymmetric case
336    returns:
337        iabsnt = True is reflection is forbidden by symmetry
338        mulp = reflection multiplicity including Fridel pairs
339        Uniq = numpy array of equivalent hkl in descending order of h,k,l
340        Phs = numpy array of corresponding phase factors (multiples of 2pi)
341       
342    NB: only works on one HKL at a time -
343        it can be made to do an array of HKLs - not any faster!
344    ''' 
345    iabsnt = False
346    H = np.array(HKL)
347    Inv = SGData['SGInv']
348    Cen = SGData['SGCen'][1:]           #skip P
349    # do cell centering extinction check
350    for cen in Cen:
351        iabsnt |= bool(int(round(np.sum(cen*H*12)))%12)
352    # get operators & expand if centrosymmetric
353    Ops = SGData['SGOps']
354    opM = np.array([op[0] for op in Ops])
355    opT = np.array([op[1] for op in Ops])
356    if Inv:
357        opM = np.concatenate((opM,-opM))
358        opT = np.concatenate((opT,-opT))
359    # generate full hkl table and phase factors; find duplicates for multiplicity
360    eqH = np.inner(opM,H)
361    HT = np.dot(opT,H)%1.
362    dup = len(ma.nonzero([ma.allclose(H,hkl,atol=0.001) for hkl in eqH])[0])
363    mulp = len(eqH)/dup
364    # generate unique reflection set (with or without Friedel pairs)
365    HPT = [tuple(hp) for hp in np.column_stack((eqH,HT))]
366    HPT.sort()
367    UPT = HPT[::dup]
368    UPT.reverse()
369    if Inv and not Friedel:
370        UPT = UPT[:len(UPT)/2]           
371    Uniq = np.split(UPT,[3,4],axis=1)
372    Phs = Uniq[1].T[0]
373    Uniq = Uniq[0]
374    # determine if extinct
375    HPT = np.split(HPT,[3,4],axis=1)[1].T[0]
376    HPT = np.reshape(HPT,(mulp,-1)).T
377    HPT = np.around(HPT-HPT[0],4)%1.
378    iabsnt |= np.any(HPT)
379    return iabsnt,mulp,Uniq
380                                   
381#def GenHKLs(H,SGData,Friedel=False):
382#    '''
383#    Generate set of equivalent reflections
384#    input:
385#        H - np.array of [h,k,l] assumed to exclude lattice centering extinctions
386#        SGData - space group data obtained from SpcGroup
387#        Friedel = True to retain Friedel pairs for centrosymmetric case
388#    returns:
389#        iabsnt = True is reflection is forbidden by symmetry
390#        mulp = reflection multiplicity including Fridel pairs
391#        Uniq = numpy array of equivalent hkl in descending order of h,k,l
392#        Phs = numpy array of corresponding phase factors (multiples of 2pi)
393#    this is not any faster than GenHKL above - oh well
394#    '''
395#    H = np.array(H,dtype=float)
396#    Inv = SGData['SGInv']
397#    Ops = SGData['SGOps']
398#    opM = np.array([op[0] for op in Ops])
399#    opT = np.array([op[1] for op in Ops])
400#    if Inv:
401#        opM = np.concatenate((opM,-opM))
402#        opT = np.concatenate((opT,-opT))
403#    # generate full hkl table and phase factors; find duplicates for multiplicity
404#    eqH = np.swapaxes(np.swapaxes(np.inner(opM,H),0,2),1,2)
405#    HeqH = zip(H,eqH)
406#    dup = np.array([sum([ma.allclose(h,hkl,atol=0.001) for hkl in HKL]) for h,HKL in HeqH])
407#    mulp = len(eqH[0])/dup
408#    # generate unique reflection set (with or without Friedel pairs)
409#    HT = np.dot(H,opT.T)%1.
410#    HPT = zip(mulp,dup,[[tuple(hp) for hp in HP] for HP in np.dstack((eqH,HT))])
411#    Uniq = []
412#    Phs = []
413#    iabsnt = []
414#    for mp,dp,hpt in HPT:
415#        hpt.sort()
416#        upt = hpt[::dp]
417#        upt.reverse()
418#        if Inv and not Friedel:
419#            upt = upt[:len(upt)/2]           
420#        uniq = np.split(upt,[3,4],axis=1)
421#        Phs.append(uniq[1].T[0])
422#        Uniq.append(uniq[0])
423#        # determine if extinct
424#        hpt = np.split(hpt,[3,4],axis=1)[1].T[0]
425#        hpt = np.reshape(hpt,(mp,-1)).T
426#        hpt = np.around(hpt-hpt[0],4)%1.
427#        iabsnt.append(np.any(hpt))
428#    return iabsnt,mulp,Uniq,Phs
429                                   
430def GetOprPtrName(key):           
431    OprPtrName = {
432        '-6643':[   2,' 1bar ', 1],'6479' :[  10,'  2z  ', 2],'-6479':[   9,'  mz  ', 3],
433        '6481' :[   7,'  my  ', 4],'-6481':[   6,'  2y  ', 5],'6641' :[   4,'  mx  ', 6],
434        '-6641':[   3,'  2x  ', 7],'6591' :[  28,' m+-0 ', 8],'-6591':[  27,' 2+-0 ', 9],
435        '6531' :[  25,' m110 ',10],'-6531':[  24,' 2110 ',11],'6537' :[  61,'  4z  ',12],
436        '-6537':[  62,' -4z  ',13],'975'  :[  68,' 3+++1',14],'6456' :[ 114,'  3z1 ',15],
437        '-489' :[  73,' 3+-- ',16],'483'  :[  78,' 3-+- ',17],'-969' :[  83,' 3--+ ',18],
438        '819'  :[  22,' m+0- ',19],'-819' :[  21,' 2+0- ',20],'2431' :[  16,' m0+- ',21],
439        '-2431':[  15,' 20+- ',22],'-657' :[  19,' m101 ',23],'657'  :[  18,' 2101 ',24],
440        '1943' :[  48,' -4x  ',25],'-1943':[  47,'  4x  ',26],'-2429':[  13,' m011 ',27],
441        '2429' :[  12,' 2011 ',28],'639'  :[  55,' -4y  ',29],'-639' :[  54,'  4y  ',30],
442        '-6484':[ 146,' 2010 ', 4],'6484' :[ 139,' m010 ', 5],'-6668':[ 145,' 2100 ', 6],
443        '6668' :[ 138,' m100 ', 7],'-6454':[ 148,' 2120 ',18],'6454' :[ 141,' m120 ',19],
444        '-6638':[ 149,' 2210 ',20],'6638' :[ 142,' m210 ',21],              #search ends here
445        '2223' :[  68,' 3+++2',39],
446        '6538' :[ 106,'  6z1 ',40],'-2169':[  83,' 3--+2',41],'2151' :[  73,' 3+--2',42],
447        '2205' :[  79,'-3-+-2',43],'-2205':[  78,' 3-+-2',44],'489'  :[  74,'-3+--1',45],
448        '801'  :[  53,'  4y1 ',46],'1945' :[  47,'  4x3 ',47],'-6585':[  62,' -4z3 ',48],
449        '6585' :[  61,'  4z3 ',49],'6584' :[ 114,'  3z2 ',50],'6666' :[ 106,'  6z5 ',51],
450        '6643' :[   1,' Iden ',52],'-801' :[  55,' -4y1 ',53],'-1945':[  48,' -4x3 ',54],
451        '-6666':[ 105,' -6z5 ',55],'-6538':[ 105,' -6z1 ',56],'-2223':[  69,'-3+++2',57],
452        '-975' :[  69,'-3+++1',58],'-6456':[ 113,' -3z1 ',59],'-483' :[  79,'-3-+-1',60],
453        '969'  :[  84,'-3--+1',61],'-6584':[ 113,' -3z2 ',62],'2169' :[  84,'-3--+2',63],
454        '-2151':[  74,'-3+--2',64],'0':[0,' ????',0]
455        }
456    return OprPtrName[key]
457
458def GetKNsym(key):
459    KNsym = {
460        '0'         :'    1   ','1'         :'   -1   ','64'        :'  2(100)','32'        :'  m(100)',
461        '97'        :'2/m(100)','16'        :'  2(010)','8'         :'  m(010)','25'        :'2/m(010)',
462        '2'         :'  2(001)','4'         :'  m(001)','7'         :'2/m(001)','134217728' :'  2(011)',
463        '67108864'  :'  m(011)','201326593' :'2/m(011)','2097152'   :'  2(0+-)','1048576'   :'  m(0+-)',
464        '3145729'   :'2/m(0+-)','8388608'   :'  2(101)','4194304'   :'  m(101)','12582913'  :'2/m(101)',
465        '524288'    :'  2(+0-)','262144'    :'  m(+0-)','796433'    :'2/m(+0-)','1024'      :'  2(110)',
466        '512'       :'  m(110)','1537'      :'2/m(110)','256'       :'  2(+-0)','128'       :'  m(+-0)',
467        '385'       :'2/m(+-0)','76'        :'mm2(100)','52'        :'mm2(010)','42'        :'mm2(001)',
468        '135266336' :'mm2(011)','69206048'  :'mm2(0+-)','8650760'   :'mm2(101)','4718600'   :'mm2(+0-)',
469        '1156'      :'mm2(110)','772'       :'mm2(+-0)','82'        :'  222   ','136314944' :'222(100)',
470        '8912912'   :'222(010)','1282'      :'222(001)','127'       :'  mmm   ','204472417' :'mmm(100)',
471        '13369369'  :'mmm(010)','1927'      :'mmm(001)','33554496'  :'  4(100)','16777280'  :' -4(100)',
472        '50331745'  :'4/m(100)','169869394' :'422(100)','84934738'  :'-42m 100','101711948' :'4mm(100)',
473        '254804095' :'4/mmm100','536870928 ':'  4(010)','268435472' :' -4(010)','805306393' :'4/m (10)',
474        '545783890' :'422(010)','272891986' :'-42m 010','541327412' :'4mm(010)','818675839' :'4/mmm010',
475        '2050'      :'  4(001)','4098'      :' -4(001)','6151'      :'4/m(001)','3410'      :'422(001)',
476        '4818'      :'-42m 001','2730'      :'4mm(001)','8191'      :'4/mmm001','8192'      :'  3(111)',
477        '8193'      :' -3(111)','2629888'   :' 32(111)','1319040'   :' 3m(111)','3940737'   :'-3m(111)',
478        '32768'     :'  3(+--)','32769'     :' -3(+--)','10519552'  :' 32(+--)','5276160'   :' 3m(+--)',
479        '15762945'  :'-3m(+--)','65536'     :'  3(-+-)','65537'     :' -3(-+-)','134808576' :' 32(-+-)',
480        '67437056'  :' 3m(-+-)','202180097' :'-3m(-+-)','131072'    :'  3(--+)','131073'    :' -3(--+)',
481        '142737664' :' 32(--+)','71434368'  :' 3m(--+)','214040961' :'-3m(--+)','237650'    :'   23   ',
482        '237695'    :'   m3   ','715894098' :'   432  ','358068946' :'  -43m  ','1073725439':'   m3m  ',
483        '68157504'  :' mm2d100','4456464'   :' mm2d010','642'       :' mm2d001','153092172' :'-4m2 100',
484        '277348404' :'-4m2 010','5418'      :'-4m2 001','1075726335':'  6/mmm ','1074414420':'-6m2 100',
485        '1075070124':'-6m2 120','1075069650':'   6mm  ','1074414890':'   622  ','1073758215':'   6/m  ',
486        '1073758212':'   -6   ','1073758210':'    6   ','1073759865':'-3m(100)','1075724673':'-3m(120)',
487        '1073758800':' 3m(100)','1075069056':' 3m(120)','1073759272':' 32(100)','1074413824':' 32(120)',
488        '1073758209':'   -3   ','1073758208':'    3   ','1074135143':'mmm(100)','1075314719':'mmm(010)',
489        '1073743751':'mmm(110)','1074004034':' mm2z100','1074790418':' mm2z010','1073742466':' mm2z110',
490        '1074004004':'mm2(100)','1074790412':'mm2(010)','1073742980':'mm2(110)','1073872964':'mm2(120)',
491        '1074266132':'mm2(210)','1073742596':'mm2(+-0)','1073872930':'222(100)','1074266122':'222(010)',
492        '1073743106':'222(110)','1073741831':'2/m(001)','1073741921':'2/m(100)','1073741849':'2/m(010)',
493        '1073743361':'2/m(110)','1074135041':'2/m(120)','1075314689':'2/m(210)','1073742209':'2/m(+-0)',
494        '1073741828':' m(001) ','1073741888':' m(100) ','1073741840':' m(010) ','1073742336':' m(110) ',
495        '1074003968':' m(120) ','1074790400':' m(210) ','1073741952':' m(+-0) ','1073741826':' 2(001) ',
496        '1073741856':' 2(100) ','1073741832':' 2(010) ','1073742848':' 2(110) ','1073872896':' 2(120) ',
497        '1074266112':' 2(210) ','1073742080':' 2(+-0) ','1073741825':'   -1   '
498        }
499    return KNsym[key]       
500
501def GetNXUPQsym(siteSym):       
502    NXUPQsym = {
503        '    1   ':(28,29,28,28),'   -1   ':( 1,29,28, 0),'  2(100)':(12,18,12,25),'  m(100)':(25,18,12,25),
504        '2/m(100)':( 1,18, 0,-1),'  2(010)':(13,17,13,24),'  m(010)':(24,17,13,24),'2/m(010)':( 1,17, 0,-1),
505        '  2(001)':(14,16,14,23),'  m(001)':(23,16,14,23),'2/m(001)':( 1,16, 0,-1),'  2(011)':(10,23,10,22),
506        '  m(011)':(22,23,10,22),'2/m(011)':( 1,23, 0,-1),'  2(0+-)':(11,24,11,21),'  m(0+-)':(21,24,11,21),
507        '2/m(0+-)':( 1,24, 0,-1),'  2(101)':( 8,21, 8,20),'  m(101)':(20,21, 8,20),'2/m(101)':( 1,21, 0,-1),
508        '  2(+0-)':( 9,22, 9,19),'  m(+0-)':(19,22, 9,19),'2/m(+0-)':( 1,22, 0,-1),'  2(110)':( 6,19, 6,18),
509        '  m(110)':(18,19, 6,18),'2/m(110)':( 1,19, 0,-1),'  2(+-0)':( 7,20, 7,17),'  m(+-0)':(17,20, 7,17),
510        '2/m(+-0)':( 1,20, 0,-1),'mm2(100)':(12,10, 0,-1),'mm2(010)':(13,10, 0,-1),'mm2(001)':(14,10, 0,-1),
511        'mm2(011)':(10,13, 0,-1),'mm2(0+-)':(11,13, 0,-1),'mm2(101)':( 8,12, 0,-1),'mm2(+0-)':( 9,12, 0,-1),
512        'mm2(110)':( 6,11, 0,-1),'mm2(+-0)':( 7,11, 0,-1),'  222   ':( 1,10, 0,-1),'222(100)':( 1,13, 0,-1),
513        '222(010)':( 1,12, 0,-1),'222(001)':( 1,11, 0,-1),'  mmm   ':( 1,10, 0,-1),'mmm(100)':( 1,13, 0,-1),
514        'mmm(010)':( 1,12, 0,-1),'mmm(001)':( 1,11, 0,-1),'  4(100)':(12, 4,12, 0),' -4(100)':( 1, 4,12, 0),
515        '4/m(100)':( 1, 4,12,-1),'422(100)':( 1, 4, 0,-1),'-42m 100':( 1, 4, 0,-1),'4mm(100)':(12, 4, 0,-1),
516        '4/mmm100':( 1, 4, 0,-1),'  4(010)':(13, 3,13, 0),' -4(010)':( 1, 3,13, 0),'4/m (10)':( 1, 3,13,-1),
517        '422(010)':( 1, 3, 0,-1),'-42m 010':( 1, 3, 0,-1),'4mm(010)':(13, 3, 0,-1),'4/mmm010':(1, 3, 0,-1,),
518        '  4(001)':(14, 2,14, 0),' -4(001)':( 1, 2,14, 0),'4/m(001)':( 1, 2,14,-1),'422(001)':( 1, 2, 0,-1),
519        '-42m 001':( 1, 2, 0,-1),'4mm(001)':(14, 2, 0,-1),'4/mmm001':( 1, 2, 0,-1),'  3(111)':( 2, 5, 2, 0),
520        ' -3(111)':( 1, 5, 2, 0),' 32(111)':( 1, 5, 0, 2),' 3m(111)':( 2, 5, 0, 2),'-3m(111)':( 1, 5, 0,-1),
521        '  3(+--)':( 5, 8, 5, 0),' -3(+--)':( 1, 8, 5, 0),' 32(+--)':( 1, 8, 0, 5),' 3m(+--)':( 5, 8, 0, 5),
522        '-3m(+--)':( 1, 8, 0,-1),'  3(-+-)':( 4, 7, 4, 0),' -3(-+-)':( 1, 7, 4, 0),' 32(-+-)':( 1, 7, 0, 4),
523        ' 3m(-+-)':( 4, 7, 0, 4),'-3m(-+-)':( 1, 7, 0,-1),'  3(--+)':( 3, 6, 3, 0),' -3(--+)':( 1, 6, 3, 0),
524        ' 32(--+)':( 1, 6, 0, 3),' 3m(--+)':( 3, 6, 0, 3),'-3m(--+)':( 1, 6, 0,-1),'   23   ':( 1, 1, 0, 0),
525        '   m3   ':( 1, 1, 0, 0),'   432  ':( 1, 1, 0, 0),'  -43m  ':( 1, 1, 0, 0),'   m3m  ':( 1, 1, 0, 0),
526        ' mm2d100':(12,13, 0,-1),' mm2d010':(13,12, 0,-1),' mm2d001':(14,11, 0,-1),'-4m2 100':( 1, 4, 0,-1),
527        '-4m2 010':( 1, 3, 0,-1),'-4m2 001':( 1, 2, 0,-1),'  6/mmm ':( 1, 9, 0,-1),'-6m2 100':( 1, 9, 0,-1),
528        '-6m2 120':( 1, 9, 0,-1),'   6mm  ':(14, 9, 0,-1),'   622  ':( 1, 9, 0,-1),'   6/m  ':( 1, 9,14,-1),
529        '   -6   ':( 1, 9,14, 0),'    6   ':(14, 9,14, 0),'-3m(100)':( 1, 9, 0,-1),'-3m(120)':( 1, 9, 0,-1),
530        ' 3m(100)':(14, 9, 0,14),' 3m(120)':(14, 9, 0,14),' 32(100)':( 1, 9, 0,14),' 32(120)':( 1, 9, 0,14),
531        '   -3   ':( 1, 9,14, 0),'    3   ':(14, 9,14, 0),'mmm(100)':( 1,14, 0,-1),'mmm(010)':( 1,15, 0,-1),
532        'mmm(110)':( 1,11, 0,-1),' mm2z100':(14,14, 0,-1),' mm2z010':(14,15, 0,-1),' mm2z110':(14,11, 0,-1),
533        'mm2(100)':(12,14, 0,-1),'mm2(010)':(13,15, 0,-1),'mm2(110)':( 6,11, 0,-1),'mm2(120)':(15,14, 0,-1),
534        'mm2(210)':(16,15, 0,-1),'mm2(+-0)':( 7,11, 0,-1),'222(100)':( 1,14, 0,-1),'222(010)':( 1,15, 0,-1),
535        '222(110)':( 1,11, 0,-1),'2/m(001)':( 1,16,14,-1),'2/m(100)':( 1,25,12,-1),'2/m(010)':( 1,28,13,-1),
536        '2/m(110)':( 1,19, 6,-1),'2/m(120)':( 1,27,15,-1),'2/m(210)':( 1,26,16,-1),'2/m(+-0)':( 1,20,17,-1),
537        ' m(001) ':(23,16,14,23),' m(100) ':(26,25,12,26),' m(010) ':(27,28,13,27),' m(110) ':(18,19, 6,18),
538        ' m(120) ':(24,27,15,24),' m(210) ':(25,26,16,25),' m(+-0) ':(17,20, 7,17),' 2(001) ':(14,16,14,23),
539        ' 2(100) ':(12,25,12,26),' 2(010) ':(13,28,13,27),' 2(110) ':( 6,19, 6,18),' 2(120) ':(15,27,15,24),
540        ' 2(210) ':(16,26,16,25),' 2(+-0) ':( 7,20, 7,17),'   -1   ':( 1,29,28, 0)
541        }
542    return NXUPQsym[siteSym]
543
544def GetCSxinel(siteSym): 
545    CSxinel = [[],                         # 0th empty - indices are Fortran style
546        [[0,0,0],[ 0.0, 0.0, 0.0]],      #  0  0  0
547        [[1,1,1],[ 1.0, 1.0, 1.0]],      #  X  X  X
548        [[1,1,1],[ 1.0, 1.0,-1.0]],      #  X  X -X
549        [[1,1,1],[ 1.0,-1.0, 1.0]],      #  X -X  X
550        [[1,1,1],[ 1.0,-1.0,-1.0]],      # -X  X  X
551        [[1,1,0],[ 1.0, 1.0, 0.0]],      #  X  X  0
552        [[1,1,0],[ 1.0,-1.0, 0.0]],      #  X -X  0
553        [[1,0,1],[ 1.0, 0.0, 1.0]],      #  X  0  X
554        [[1,0,1],[ 1.0, 0.0,-1.0]],      #  X  0 -X
555        [[0,1,1],[ 0.0, 1.0, 1.0]],      #  0  Y  Y
556        [[0,1,1],[ 0.0, 1.0,-1.0]],      #  0  Y -Y
557        [[1,0,0],[ 1.0, 0.0, 0.0]],      #  X  0  0
558        [[0,1,0],[ 0.0, 1.0, 0.0]],      #  0  Y  0
559        [[0,0,1],[ 0.0, 0.0, 1.0]],      #  0  0  Z
560        [[1,1,0],[ 1.0, 2.0, 0.0]],      #  X 2X  0
561        [[1,1,0],[ 2.0, 1.0, 0.0]],      # 2X  X  0
562        [[1,1,2],[ 1.0, 1.0, 1.0]],      #  X  X  Z
563        [[1,1,2],[ 1.0,-1.0, 1.0]],      #  X -X  Z
564        [[1,2,1],[ 1.0, 1.0, 1.0]],      #  X  Y  X
565        [[1,2,1],[ 1.0, 1.0,-1.0]],      #  X  Y -X
566        [[1,2,2],[ 1.0, 1.0, 1.0]],      #  X  Y  Y
567        [[1,2,2],[ 1.0, 1.0,-1.0]],      #  X  Y -Y
568        [[1,2,0],[ 1.0, 1.0, 0.0]],      #  X  Y  0
569        [[1,0,2],[ 1.0, 0.0, 1.0]],      #  X  0  Z
570        [[0,1,2],[ 0.0, 1.0, 1.0]],      #  0  Y  Z
571        [[1,1,2],[ 1.0, 2.0, 1.0]],      #  X 2X  Z
572        [[1,1,2],[ 2.0, 1.0, 1.0]],      # 2X  X  Z
573        [[1,2,3],[ 1.0, 1.0, 1.0]],      #  X  Y  Z
574        ]
575    indx = GetNXUPQsym(siteSym)
576    return CSxinel[indx[0]]
577   
578def GetCSuinel(siteSym):
579    # returns Uij terms, multipliers, GUI flags & Uiso2Uij multipliers
580    CSuinel = [[],                                             # 0th empty - indices are Fortran style
581        [[1,1,1,0,0,0],[ 1.0, 1.0, 1.0, 0.0, 0.0, 0.0],[1,0,0,0,0,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  A  A  0  0  0
582        [[1,1,2,0,0,0],[ 1.0, 1.0, 1.0, 0.0, 0.0, 0.0],[1,0,1,0,0,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  A  C  0  0  0
583        [[1,2,1,0,0,0],[ 1.0, 1.0, 1.0, 0.0, 0.0, 0.0],[1,1,0,0,0,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  A  0  0  0
584        [[1,2,2,0,0,0],[ 1.0, 1.0, 1.0, 0.0, 0.0, 0.0],[1,1,0,0,0,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  B  0  0  0
585        [[1,1,1,2,2,2],[ 1.0, 1.0, 1.0, 1.0, 1.0, 1.0],[1,0,0,1,0,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  A  A  D  D  D
586        [[1,1,1,2,2,2],[ 1.0, 1.0, 1.0, 1.0,-1.0,-1.0],[1,0,0,1,0,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  A  A  D -D -D
587        [[1,1,1,2,2,2],[ 1.0, 1.0, 1.0, 1.0,-1.0, 1.0],[1,0,0,1,0,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  A  A  D -D  D
588        [[1,1,1,2,2,2],[ 1.0, 1.0, 1.0, 1.0, 1.0,-1.0],[1,0,0,1,0,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  A  A  D  D -D
589        [[1,1,2,1,0,0],[ 1.0, 1.0, 1.0, 0.5, 0.0, 0.0],[1,0,1,0,0,0],[1.0,1.0,1.0,0.5,0.0,0.0]],    #  A  A  C A/2 0  0
590        [[1,2,3,0,0,0],[ 1.0, 1.0, 1.0, 0.0, 0.0, 0.0],[1,1,1,0,0,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  C  0  0  0
591        [[1,1,2,3,0,0],[ 1.0, 1.0, 1.0, 1.0, 0.0, 0.0],[1,0,1,1,0,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  A  C  D  0  0
592        [[1,2,1,0,3,0],[ 1.0, 1.0, 1.0, 0.0, 1.0, 0.0],[1,1,0,0,1,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  A  0  E  0
593        [[1,2,2,0,0,3],[ 1.0, 1.0, 1.0, 0.0, 0.0, 1.0],[1,1,0,0,0,1],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  B  0  0  F
594        [[1,2,3,2,0,0],[ 1.0, 1.0, 1.0, 0.5, 0.0, 0.0],[1,1,1,0,0,0],[1.0,1.0,1.0,0.0,0.5,0.0]],    #  A  B  C B/2 0  0
595        [[1,2,3,1,0,0],[ 1.0, 1.0, 1.0, 0.5, 0.0, 0.0],[1,1,1,0,0,0],[1.0,1.0,1.0,0.0,0.5,0.0]],    #  A  B  C A/2 0  0
596        [[1,2,3,4,0,0],[ 1.0, 1.0, 1.0, 1.0, 0.0, 0.0],[1,1,1,1,0,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  C  D  0  0
597        [[1,2,3,0,4,0],[ 1.0, 1.0, 1.0, 0.0, 1.0, 0.0],[1,1,1,0,1,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  C  0  E  0
598        [[1,2,3,0,0,4],[ 1.0, 1.0, 1.0, 0.0, 0.0, 1.0],[1,1,1,0,0,1],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  C  0  0  F
599        [[1,1,2,3,4,4],[ 1.0, 1.0, 1.0, 1.0, 1.0,-1.0],[1,0,1,1,1,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  A  C  D  E -E
600        [[1,1,2,3,4,4],[ 1.0, 1.0, 1.0, 1.0, 1.0, 1.0],[1,0,1,1,1,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  A  C  D  E  E
601        [[1,2,1,3,4,3],[ 1.0, 1.0, 1.0, 1.0, 1.0,-1.0],[1,1,0,1,1,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  A  D  E -D
602        [[1,2,1,3,4,3],[ 1.0, 1.0, 1.0, 1.0, 1.0, 1.0],[1,1,0,1,1,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  A  D  E  D
603        [[1,2,2,3,3,4],[ 1.0, 1.0, 1.0, 1.0,-1.0, 1.0],[1,1,0,1,0,1],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  B  D -D  F
604        [[1,2,2,3,3,4],[ 1.0, 1.0, 1.0, 1.0, 1.0, 1.0],[1,1,0,1,0,1],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  B  D  D  F
605        [[1,2,3,2,4,4],[ 1.0, 1.0, 1.0, 0.5, 0.5, 1.0],[1,1,1,0,0,1],[1.0,1.0,1.0,0.5,0.0,0.0]],    #  A  B  C B/2 F/2 F
606        [[1,2,3,1,0,4],[ 1.0, 1.0, 1.0, 0.5, 0.0, 1.0],[1,1,1,0,0,1],[1.0,1.0,1.0,0.5,0.0,0.0]],    #  A  B  C A/2  0  F
607        [[1,2,3,2,4,0],[ 1.0, 1.0, 1.0, 0.5, 1.0, 0.0],[1,1,1,0,1,0],[1.0,1.0,1.0,0.5,0.0,0.0]],    #  A  B  C B/2  E  0
608        [[1,2,3,1,4,4],[ 1.0, 1.0, 1.0, 0.5, 1.0, 0.5],[1,1,1,0,1,0],[1.0,1.0,1.0,0.5,0.0,0.0]],    #  A  B  C A/2  E E/2
609        [[1,2,3,4,5,6],[ 1.0, 1.0, 1.0, 1.0, 1.0, 1.0],[1,1,1,1,1,1],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  C  D  E   F
610        ]
611    indx = GetNXUPQsym(siteSym)
612    return CSuinel[indx[1]]
613   
614def MustrainNames(SGData):
615    laue = SGData['SGLaue']
616    uniq = SGData['SGUniq']
617    if laue in ['m3','m3m']:
618        return ['S400','S220']
619    elif laue in ['6/m','6/mmm','3m1']:
620        return ['S400','S004','S202']
621    elif laue in ['31m','3']:
622        return ['S400','S004','S202','S211']
623    elif laue in ['3R','3mR']:
624        return ['S400','S220','S310','S211']
625    elif laue in ['4/m','4/mmm']:
626        return ['S400','S004','S220','S022']
627    elif laue in ['mmm']:
628        return ['S400','S040','S004','S220','S202','S022']
629    elif laue in ['2/m']:
630        SHKL = ['S400','S040','S004','S220','S202','S022']
631        if uniq == 'a':
632            SHKL += ['S013','S031','S211']
633        elif uniq == 'b':
634            SHKL += ['S301','S103','S121']
635        elif uniq == 'c':
636            SHKL += ['S130','S310','S112']
637        return SHKL
638    else:
639        SHKL = ['S400','S040','S004','S220','S202','S022']
640        SHKL += ['S310','S103','S031','S130','S301','S013']
641        SHKL += ['S211','S121','S112']
642        return SHKL
643
644def MustrainCoeff(HKL,SGData):
645    #NB: order of terms is the same as returned by MustrainNames
646    laue = SGData['SGLaue']
647    uniq = SGData['SGUniq']
648    h,k,l = HKL
649    Strm = []
650    if laue in ['m3','m3m']:
651        Strm.append(h**4+k**4+l**4)
652        Strm.append(3.0*((h*k)**2+(h*l)**2+(k*l)**2))
653    elif laue in ['6/m','6/mmm','3m1']:
654        Strm.append(h**4+k**4+2.0*k*h**3+2.0*h*k**3+3.0*(h*k)**2)
655        Strm.append(l**4)
656        Strm.append(3.0*((h*l)**2+(k*l)**2+h*k*l**2))
657    elif laue in ['31m','3']:
658        Strm.append(h**4+k**4+2.0*k*h**3+2.0*h*k**3+3.0*(h*k)**2)
659        Strm.append(l**4)
660        Strm.append(3.0*((h*l)**2+(k*l)**2+h*k*l**2))
661        Strm.append(4.0*h*k*l*(h+k))
662    elif laue in ['3R','3mR']:
663        Strm.append(h**4+k**4+l**4)
664        Strm.append(3.0*((h*k)**2+(h*l)**2+(k*l)**2))
665        Strm.append(2.0*(h*l**3+l*k**3+k*h**3)+2.0*(l*h**3+k*l**3+l*k**3))
666        Strm.append(4.0*(k*l*h**2+h*l*k**2+h*k*l**2))
667    elif laue in ['4/m','4/mmm']:
668        Strm.append(h**4+k**4)
669        Strm.append(l**4)
670        Strm.append(3.0*(h*k)**2)
671        Strm.append(3.0*((h*l)**2+(k*l)**2))
672    elif laue in ['mmm']:
673        Strm.append(h**4)
674        Strm.append(k**4)
675        Strm.append(l**4)
676        Strm.append(3.0*(h*k)**2)
677        Strm.append(3.0*(h*l)**2)
678        Strm.append(3.0*(k*l)**2)
679    elif laue in ['2/m']:
680        Strm.append(h**4)
681        Strm.append(k**4)
682        Strm.append(l**4)
683        Strm.append(3.0*(h*k)**2)
684        Strm.append(3.0*(h*l)**2)
685        Strm.append(3.0*(k*l)**2)
686        if uniq == 'a':
687            Strm.append(2.0*k*l**3)
688            Strm.append(2.0*l*k**3)
689            Strm.append(4.0*k*l*h**2)
690        elif uniq == 'b':
691            Strm.append(2.0*l*h**3)
692            Strm.append(2.0*h*l**3)
693            Strm.append(4.0*h*l*k**2)
694        elif uniq == 'c':
695            Strm.append(2.0*h*k**3)
696            Strm.append(2.0*k*h**3)
697            Strm.append(4.0*h*k*l**2)
698    else:
699        Strm.append(h**4)
700        Strm.append(k**4)
701        Strm.append(l**4)
702        Strm.append(3.0*(h*k)**2)
703        Strm.append(3.0*(h*l)**2)
704        Strm.append(3.0*(k*l)**2)
705        Strm.append(2.0*k*h**3)
706        Strm.append(2.0*h*l**3)
707        Strm.append(2.0*l*k**3)
708        Strm.append(2.0*h*k**3)
709        Strm.append(2.0*l*h**3)
710        Strm.append(2.0*k*l**3)
711        Strm.append(4.0*k*l*h**2)
712        Strm.append(4.0*h*l*k**2)
713        Strm.append(4.0*k*h*l**2)
714    return Strm
715   
716def Muiso2Shkl(muiso,SGData,cell):
717    #this is to convert isotropic mustrain to generalized Shkls - doesn't work just now
718    import GSASIIlattice as G2lat
719    from scipy.optimize import fmin
720    A = G2lat.cell2AB(cell)[0]
721    def minMus(Shkl,H,muiso,SGData,A):
722        U = np.inner(A.T,H)
723        S = np.array(MustrainCoeff(H.T,SGData))
724        sum = np.sqrt(np.sum(np.multiply(S,Shkl)))
725        return abs(muiso-sum*H)
726    laue = SGData['SGLaue']
727    if laue in ['m3','m3m']:
728        H = [[1,0,0],[1,1,0]]
729        S0 = [0.01,0.01]
730    elif laue in ['6/m','6/mmm','3m1']:
731        H = [[1,0,0],[0,0,1],[1,0,1]]
732        S0 = [0.01,0.01,0.01]
733    elif laue in ['31m','3']:
734        H = [[1,0,0],[0,0,1],[1,0,1],[1,1,1]]
735        S0 = [0.01,0.01,0.01,0.01]
736    elif laue in ['3R','3mR']:
737        H = [[1,0,0],[1,1,0],[1,0,1],[1,1,1]]
738        S0 = [0.01,0.01,0.01,0.01]
739    elif laue in ['4/m','4/mmm']:
740        H = [[1,0,0],[0,0,1],[1,1,0],[1,0,1]]
741        S0 = [0.01,0.01,0.01,0.01]
742    elif laue in ['mmm']:
743        H = [[1,0,0],[0,1,0],[0,0,1],[1,1,0],[1,0,1],[0,1,1]]
744        S0 = [0.01,0.01,0.01,0.01,0.01,0.01]
745    elif laue in ['2/m']:
746        H = [[1,0,0],[0,1,0],[0,0,1],[1,1,0],[1,0,1],[0,1,1]]
747        if uniq == 'a':
748            H.append([0,1,-1])
749            H.append([0,-2,1])
750        elif uniq == 'b':
751            H.append([1,0,-1])
752            H.append([-2,0,1])
753        elif uniq == 'c':
754            H.append([1,-1,0])
755            H.append([-2,1,0])
756        H.append([1,1,1])
757        S0 = [9*[0.01,]]
758    else:
759        H = [[1,0,0],[0,1,0],[0,0,1],[1,1,0],[1,0,1],[0,1,1],
760            [-1,1,0],[1,0,-1],[0,-1,1],[1,-2,0],[-2,0,1],[0,1,-2],
761            [1,-1,1],[-1, 1, 1],[1,-1,1]]
762        S0 = [15*[0.01,]]
763    H = np.array(H)
764    S0 = np.array(S0)
765    return fmin(minMus,S0,(H,muiso,SGData,A))
766       
767def SytSym(XYZ,SGData):
768    '''
769    Generates the number of equivalent positions and a site symmetry code for a specified coordinate and space group
770    input: 
771       XYZ: an array, tuple or list containing 3 elements: x, y & z
772       SGData: from SpcGroup
773    Returns a two element tuple:
774       The 1st element is a code for the site symmetry (see GetKNsym)
775       The 2nd element is the site multiplicity
776    '''
777    def PackRot(SGOps):
778        IRT = []
779        for ops in SGOps:
780            M = ops[0]
781            irt = 0
782            for j in range(2,-1,-1):
783                for k in range(2,-1,-1):
784                    irt *= 3
785                    irt += M[k][j]
786            IRT.append(int(irt))
787        return IRT
788       
789    SymName = ''
790    Mult = 1
791    Isym = 0
792    if SGData['SGLaue'] in ['3','3m1','31m','6/m','6/mmm']:
793        Isym = 1073741824
794    Jdup = 0
795    Xeqv = GenAtom(XYZ,SGData,True)
796    IRT = PackRot(SGData['SGOps'])
797    L = -1
798    for ic,cen in enumerate(SGData['SGCen']):
799        for invers in range(int(SGData['SGInv']+1)):
800            for io,ops in enumerate(SGData['SGOps']):
801                irtx = (1-2*invers)*IRT[io]
802                L += 1
803                if not Xeqv[L][1]:
804                    Jdup += 1
805                    jx = GetOprPtrName(str(irtx))
806                    if jx[2] < 39:
807                        Isym += 2**(jx[2]-1)
808    if Isym == 1073741824: Isym = 0
809    Mult = len(SGData['SGOps'])*len(SGData['SGCen'])*(int(SGData['SGInv'])+1)/Jdup
810         
811    return GetKNsym(str(Isym)),Mult
812   
813def ElemPosition(SGData):
814    ''' Under development
815    Object here is to return a list of symmetry element types and locations suitable
816    for say drawing them.
817    So far I have the element type... getting all possible locations without lookup may be impossible!
818    '''
819    SymElements = []
820    Inv = SGData['SGInv']
821    Cen = SGData['SGCen']
822    eleSym = {-3:['','-1'],-2:['',-6],-1:['2','-4'],0:['3','-3'],1:['4','m'],2:['6',''],3:['1','']}
823    # get operators & expand if centrosymmetric
824    Ops = SGData['SGOps']
825    opM = np.array([op[0].T for op in Ops])
826    opT = np.array([op[1] for op in Ops])
827    if Inv:
828        opM = np.concatenate((opM,-opM))
829        opT = np.concatenate((opT,-opT))
830    opMT = zip(opM,opT)
831    for M,T in opMT[1:]:        #skip I
832        Dt = int(nl.det(M))
833        Tr = int(np.trace(M))
834        Dt = -(Dt-1)/2
835        Es = eleSym[Tr][Dt]
836        if Dt:              #rotation-inversion
837            I = np.eye(3)
838            if Tr == 1:     #mirrors/glides
839                if np.any(T):       #glide
840                    M2 = np.inner(M,M)
841                    MT = np.inner(M,T)+T
842                    print 'glide',Es,MT
843                    print M2
844                else:               #mirror
845                    print 'mirror',Es,T
846                    print I-M
847                X = [-1,-1,-1]
848            elif Tr == -3:  # pure inversion
849                X = np.inner(nl.inv(I-M),T)
850                print 'inversion',Es,X
851            else:           #other rotation-inversion
852                M2 = np.inner(M,M)
853                MT = np.inner(M,T)+T
854                print 'rot-inv',Es,MT
855                print M2
856                X = [-1,-1,-1]
857               
858           
859           
860        else:               #rotations
861            print 'rotation',Es
862            X = [-1,-1,-1]
863        #SymElements.append([Es,X])
864       
865    return #SymElements
866   
867def ApplyStringOps(A,SGData,X,Uij=[]):
868    SGOps = SGData['SGOps']
869    SGCen = SGData['SGCen']
870    Ax = A.split('+')
871    Ax[0] = int(Ax[0])
872    iC = 0
873    if Ax[0] < 0:
874        iC = 1
875    Ax[0] = abs(Ax[0])
876    nA = Ax[0]%100-1
877    cA = Ax[0]/100
878    Cen = SGCen[cA]
879    M,T = SGOps[nA]
880    if len(Ax)>1:
881        cellA = Ax[1].split(',')
882        cellA = np.array([int(a) for a in cellA])
883    else:
884        cellA = np.zeros(3)
885    newX = (1-2*iC)*(Cen+np.inner(M,X)+T)+cellA
886    if len(Uij):
887        U = Uij2U(Uij)
888        U = np.inner(M,np.inner(U,M).T)
889        newUij = U2Uij(U)
890        return [newX,newUij]
891    else:
892        return newX
893       
894def StringOpsProd(A,B,SGData):
895    ''' Find A*B where A & B are in strings '-' + '100*c+n' + '+ijk'
896    where '-' indicates inversion, c(>0) is the cell centering operator,
897    n is operator number from SgOps and ijk are unit cell translations (each may be <0).
898    Should return resultant string - C.
899        SGData - dictionary using entries:
900             'SGCen': cell centering vectors [0,0,0] at least
901             'SGOps': symmetry operations as [M,T] so that M*x+T = x'
902    '''
903    SGOps = SGData['SGOps']
904    SGCen = SGData['SGCen']
905    #1st split out the cell translation part & work on the operator parts
906    Ax = A.split('+'); Bx = B.split('+')
907    Ax[0] = int(Ax[0]); Bx[0] = int(Bx[0])
908    iC = 0
909    if Ax[0]*Bx[0] < 0:
910        iC = 1
911    Ax[0] = abs(Ax[0]); Bx[0] = abs(Bx[0])
912    nA = Ax[0]%100-1;  nB = Bx[0]%100-1
913    cA = Ax[0]/100;  cB = Bx[0]/100
914    Cen = (SGCen[cA]+SGCen[cB])%1.0
915    cC = np.nonzero([np.allclose(C,Cen) for C in SGCen])[0][0]
916    Ma,Ta = SGOps[nA]; Mb,Tb = SGOps[nB]
917    Mc = np.inner(Ma,Mb.T)
918#    print Ma,Mb,Mc
919    Tc = (np.add(np.inner(Mb,Ta)+1.,Tb))%1.0
920#    print Ta,Tb,Tc
921#    print [np.allclose(M,Mc)&np.allclose(T,Tc) for M,T in SGOps]
922    nC = np.nonzero([np.allclose(M,Mc)&np.allclose(T,Tc) for M,T in SGOps])[0][0]
923    #now the cell translation part
924    if len(Ax)>1:
925        cellA = Ax[1].split(',')
926        cellA = [int(a) for a in cellA]
927    else:
928        cellA = [0,0,0]
929    if len(Bx)>1:
930        cellB = Bx[1].split(',')
931        cellB = [int(b) for b in cellB]
932    else:
933        cellB = [0,0,0]
934    cellC = np.add(cellA,cellB)
935    C = str(((nC+1)+(100*cC))*(1-2*iC))+'+'+ \
936        str(int(cellC[0]))+','+str(int(cellC[1]))+','+str(int(cellC[2]))
937    return C
938           
939def U2Uij(U):
940    #returns the UIJ vector U11,U22,U33,U12,U13,U23 from tensor U
941    return [U[0][0],U[1][1],U[2][2],U[0][1],U[0][2],U[1][2]]
942   
943def Uij2U(Uij):
944    #returns the thermal motion tensor U from Uij as numpy array
945    return np.array([[Uij[0],Uij[3],Uij[4]],[Uij[3],Uij[1],Uij[5]],[Uij[4],Uij[5],Uij[2]]])
946   
947'''A list of space groups as ordered and named in the pre-2002 International
948Tables Volume A, except that spaces are used following the GSAS convention to
949separate the different crystallographic directions.
950Note that the symmetry codes here will recognize many non-standard space group
951symbols with different settings.
952'''
953spglist = {
954    'triclinic' : ('P 1','P -1',), # 1-2
955    'monoclinic': ('P 2','P 21','C 2','P m','P c','C m','C c','P 2/m','P 21/m',
956                   'C 2/m','P 2/c','P 21/c','C 2/c',), #3-15
957    'orthorhombic': ('P 2 2 2','P 2 2 21','P 21 21 2','P 21 21 21','C 2 2 21',
958                     'C 2 2 2','F 2 2 2','I 2 2 2','I 21 21 21',
959                     'P m m 2','P m c 21','P c c 2','P m a 2','P c a 21',
960                     'P n c 2','P m n 21','P b a 2','P n a 21','P n n 2',
961                     'C m m 2','C m c 21','C c c 2','A m m 2','A b m 2',
962                     'A m a 2','A b a 2','F m m 2','F d d 2','I m m 2',
963                     'I b a 2','I m a 2','P m m m','P n n n','P c c m',
964                     'P b a n','P m m a','P n n a','P m n a','P c c a',
965                     'P b a m','P c c n','P b c m','P n n m','P m m n',
966                     'P b c n','P b c a','P n m a','C m c m','C m c a',
967                     'C m m m','C c c m','C m m a','C c c a','F m m m',
968                     'F d d d','I m m m','I b a m','I b c a','I m m a',), #16-74
969    'tetragonal': ('P 4','P 41','P 42','P 43','I 4','I 41','P -4','I -4',
970                   'P 4/m','P 42/m','P 4/n','P 42/n','I 4/m','I 41/a',
971                   'P 4 2 2','P 4 21 2','P 41 2 2','P 41 21 2','P 42 2 2',
972                   'P 42 21 2','P 43 2 2','P 43 21 2','I 4 2 2','I 41 2 2',
973                   'P 4 m m','P 4 b m','P 42 c m','P 42 n m','P 4 c c',
974                   'P 4 n c','P 42 m c','P 42 b c','I 4 m m','I 4 c m',
975                   'I 41 m d','I 41 c d','P -4 2 m','P -4 2 c','P -4 21 m',
976                   'P -4 21 c','P -4 m 2','P -4 c 2','P -4 b 2','P -4 n 2',
977                   'I -4 m 2','I -4 c 2','I -4 2 m','I -4 2 d','P 4/m m m',
978                   'P 4/m c c','P 4/n b m','P 4/n n c','P 4/m b m','P 4/m n c',
979                   'P 4/n m m','P 4/n c c','P 42/m m c','P 42/m c m',
980                   'P 42/n b c','P 42/n n m','P 42/m b c','P 42/m n m',
981                   'P 42/n m c','P 42/n c m','I 4/m m m','I 4/m c m',
982                   'I 41/a m d','I 41/a c d',), # 75-142
983    'trigonal': ('P 3','P 31','P 32','R 3','P -3','R -3','P 3 1 2','P 3 2 1',
984                 'P 31 1 2','P 31 2 1','P 32 1 2','P 32 2 1','R 3 2', 'P 3 m 1',
985                 'P 3 1 m','P 3 c 1','P 3 1 c','R 3 m','R 3 c','P -3 1 m',
986                 'P -3 1 c','P -3 m 1','P -3 c 1','R -3 m','R -3 c',), #143-167
987    'hexagonal': ('P 6','P 61','P 65','P 62','P 64','P 63','P -6','P 6/m',
988                  'P 63/m','P 6 2 2','P 61 2 2','P 65 2 2','P 62 2 2',
989                  'P 64 2 2','P 63 2 2','P 6 m m','P 6 c c','P 63 c m',
990                  'P 63 m c','P -6 m 2','P -6 c 2','P -6 2 m','P -6 2 c',
991                  'P 6/m m m','P 6/m c c','P 63/m c m','P 63/m m c',), #144-194
992    'cubic': ('P 2 3','F 2 3','I 2 3','P 21 3','I 21 3','P m -3','P n -3',
993              'F m -3','F d -3','I m -3','P a -3','I a -3','P 4 3 2','P 42 3 2',
994              'F 4 3 2','F 41 3 2','I 4 3 2','P 43 3 2','P 41 3 2','I 41 3 2',
995              'P -4 3 m','F -4 3 m','I -4 3 m','P -4 3 n','F -4 3 c','I -4 3 d',
996              'P m -3 m','P n -3 n','P m -3 n','P n -3 m','F m -3 m','F m -3 c',
997              'F d -3 m','F d -3 c','I m -3 m','I a -3 d',), #195-230
998}
999'A few non-standard space groups for test use'
1000nonstandard_sglist = ('P 21 1 1','P 1 21 1','P 1 1 21','R 3 r','R 3 2 h', 
1001                      'R -3 r', 'R 3 2 r','R 3 m h', 'R 3 m r',
1002                      'R 3 c r','R -3 c r','R -3 m r',),
1003'''A list of orthorhombic space groups that were renamed in the 2002 Volume A,
1004along with the pre-2002 name. The e designates a double glide-plane'''
1005sgequiv_2002_orthorhombic= (('A e m 2', 'A b m 2',),
1006                            ('A e a 2', 'A b a 2',),
1007                            ('C m c e', 'C m c a',),
1008                            ('C m m e', 'C m m a',),
1009                            ('C c c e', 'C c c a'),)
1010'''Use the space groups types in this order to list the symbols in the
1011order they are listed in the International Tables, vol. A'''
1012symtypelist = ('triclinic', 'monoclinic', 'orthorhombic', 'tetragonal', 
1013               'trigonal', 'hexagonal', 'cubic')
1014
1015# self-test materials follow. Requires files in directory testinp
1016def test0():
1017    '''test #0: exercise MoveToUnitCell'''
1018    msg = "MoveToUnitCell failed"
1019    assert (MoveToUnitCell([1,2,3]) == [0,0,0]).all, msg
1020    assert (MoveToUnitCell([2,-1,-2]) == [0,0,0]).all, msg
1021    assert abs(MoveToUnitCell(np.array([-.1]))[0]-0.9) < 1e-6, msg
1022    assert abs(MoveToUnitCell(np.array([.1]))[0]-0.1) < 1e-6, msg
1023
1024def test1():
1025    ''' test #1: SpcGroup and SGPrint against previous results'''
1026    testdir = ospath.join(ospath.split(ospath.abspath( __file__ ))[0],'testinp')
1027    if ospath.exists(testdir):
1028        if testdir not in sys.path: sys.path.insert(0,testdir)
1029    import spctestinp
1030    def CompareSpcGroup(spc, referr, refdict, reflist): 
1031        'Compare output from GSASIIspc.SpcGroup with results from a previous run'
1032        # if an error is reported, the dictionary can be ignored
1033        msg0 = "CompareSpcGroup failed on space group %s" % spc
1034        result = SpcGroup(spc)
1035        if result[0] == referr and referr > 0: return True
1036        keys = result[1].keys()
1037        #print result[1]['SpGrp']
1038        msg = msg0 + " in list lengths"
1039        assert len(keys) == len(refdict.keys()), msg
1040        for key in keys:
1041            if key == 'SGOps' or  key == 'SGCen':
1042                msg = msg0 + (" in key %s length" % key)
1043                assert len(refdict[key]) == len(result[1][key]), msg
1044                for i in range(len(refdict[key])):
1045                    msg = msg0 + (" in key %s level 0" % key)
1046                    assert np.allclose(result[1][key][i][0],refdict[key][i][0]), msg
1047                    msg = msg0 + (" in key %s level 1" % key)
1048                    assert np.allclose(result[1][key][i][1],refdict[key][i][1]), msg
1049            else:
1050                msg = msg0 + (" in key %s" % key)
1051                assert result[1][key] == refdict[key], msg
1052        msg = msg0 + (" in key %s reflist" % key)
1053        #for (l1,l2) in zip(reflist, SGPrint(result[1])):
1054        #    assert l2.replace('\t','').replace(' ','') == l1.replace(' ',''), 'SGPrint ' +msg
1055        assert reflist == SGPrint(result[1]), 'SGPrint ' +msg
1056    for spc in spctestinp.SGdat:
1057        CompareSpcGroup(spc, 0, spctestinp.SGdat[spc], spctestinp.SGlist[spc] )
1058
1059def test2():
1060    ''' test #2: SpcGroup against cctbx (sgtbx) computations'''
1061    testdir = ospath.join(ospath.split(ospath.abspath( __file__ ))[0],'testinp')
1062    if ospath.exists(testdir):
1063        if testdir not in sys.path: sys.path.insert(0,testdir)
1064    import sgtbxtestinp
1065    def CompareWcctbx(spcname, cctbx_in, debug=0):
1066        'Compare output from GSASIIspc.SpcGroup with results from cctbx.sgtbx'
1067        cctbx = cctbx_in[:] # make copy so we don't delete from the original
1068        spc = (SpcGroup(spcname))[1]
1069        if debug: print spc['SpGrp']
1070        if debug: print spc['SGCen']
1071        latticetype = spcname.strip().upper()[0]
1072        # lattice type of R implies Hexagonal centering", fix the rhombohedral settings
1073        if latticetype == "R" and len(spc['SGCen']) == 1: latticetype = 'P'
1074        assert latticetype == spc['SGLatt'], "Failed: %s does not match Lattice: %s" % (spcname, spc['SGLatt'])
1075        onebar = [1]
1076        if spc['SGInv']: onebar.append(-1)
1077        for (op,off) in spc['SGOps']:
1078            for inv in onebar:
1079                for cen in spc['SGCen']:
1080                    noff = off + cen
1081                    noff = MoveToUnitCell(noff)
1082                    mult = tuple((op*inv).ravel().tolist())
1083                    if debug: print "\n%s: %s + %s" % (spcname,mult,noff)
1084                    for refop in cctbx:
1085                        if debug: print refop
1086                        # check the transform
1087                        if refop[:9] != mult: continue
1088                        if debug: print "mult match"
1089                        # check the translation
1090                        reftrans = list(refop[-3:])
1091                        reftrans = MoveToUnitCell(reftrans)
1092                        if all(abs(noff - reftrans) < 1.e-5):
1093                            cctbx.remove(refop)
1094                            break
1095                    else:
1096                        assert False, "failed on %s:\n\t %s + %s" % (spcname,mult,noff)
1097    for key in sgtbxtestinp.sgtbx:
1098        CompareWcctbx(key, sgtbxtestinp.sgtbx[key])
1099
1100def test3(): 
1101    ''' test #3: exercise SytSym (includes GetOprPtrName, GenAtom, GetKNsym)
1102     for selected space groups against info in IT Volume A '''
1103    def ExerciseSiteSym (spc, crdlist):
1104        'compare site symmetries and multiplicities for a specified space group'
1105        msg = "failed on site sym test for %s" % spc
1106        (E,S) = SpcGroup(spc)
1107        assert not E, msg
1108        for t in crdlist:
1109            symb, m = SytSym(t[0],S)
1110            if symb.strip() != t[2].strip() or m != t[1]:
1111                print spc,t[0],m,symb
1112            assert m == t[1]
1113            #assert symb.strip() == t[2].strip()
1114
1115    ExerciseSiteSym('p 1',[
1116            ((0.13,0.22,0.31),1,'1'),
1117            ((0,0,0),1,'1'),
1118            ])
1119    ExerciseSiteSym('p -1',[
1120            ((0.13,0.22,0.31),2,'1'),
1121            ((0,0.5,0),1,'-1'),
1122            ])
1123    ExerciseSiteSym('C 2/c',[
1124            ((0.13,0.22,0.31),8,'1'),
1125            ((0.0,.31,0.25),4,'2(010)'),
1126            ((0.25,.25,0.5),4,'-1'),
1127            ((0,0.5,0),4,'-1'),
1128            ])
1129    ExerciseSiteSym('p 2 2 2',[
1130            ((0.13,0.22,0.31),4,'1'),
1131            ((0,0.5,.31),2,'2(001)'),
1132            ((0.5,.31,0.5),2,'2(010)'),
1133            ((.11,0,0),2,'2(100)'),
1134            ((0,0.5,0),1,'222'),
1135            ])
1136    ExerciseSiteSym('p 4/n',[
1137            ((0.13,0.22,0.31),8,'1'),
1138            ((0.25,0.75,.31),4,'2(001)'),
1139            ((0.5,0.5,0.5),4,'-1'),
1140            ((0,0.5,0),4,'-1'),
1141            ((0.25,0.25,.31),2,'4(001)'),
1142            ((0.25,.75,0.5),2,'-4(001)'),
1143            ((0.25,.75,0.0),2,'-4(001)'),
1144            ])
1145    ExerciseSiteSym('p 31 2 1',[
1146            ((0.13,0.22,0.31),6,'1'),
1147            ((0.13,0.0,0.833333333),3,'2(100)'),
1148            ((0.13,0.13,0.),3,'2(110)'),
1149            ])
1150    ExerciseSiteSym('R 3 c',[
1151            ((0.13,0.22,0.31),18,'1'),
1152            ((0.0,0.0,0.31),6,'3'),
1153            ])
1154    ExerciseSiteSym('R 3 c R',[
1155            ((0.13,0.22,0.31),6,'1'),
1156            ((0.31,0.31,0.31),2,'3(111)'),
1157            ])
1158    ExerciseSiteSym('P 63 m c',[
1159            ((0.13,0.22,0.31),12,'1'),
1160            ((0.11,0.22,0.31),6,'m(100)'),
1161            ((0.333333,0.6666667,0.31),2,'3m(100)'),
1162            ((0,0,0.31),2,'3m(100)'),
1163            ])
1164    ExerciseSiteSym('I a -3',[
1165            ((0.13,0.22,0.31),48,'1'),
1166            ((0.11,0,0.25),24,'2(100)'),
1167            ((0.11,0.11,0.11),16,'3(111)'),
1168            ((0,0,0),8,'-3(111)'),
1169            ])
1170
1171if __name__ == '__main__':
1172    test0()
1173    test1()
1174    test2()
1175    test3()
1176    print "OK"
Note: See TracBrowser for help on using the repository browser.