source: trunk/GSASIIspc.py @ 1951

Last change on this file since 1951 was 1951, checked in by vondreele, 8 years ago

implement drawing of incommensurate structures - live variation wrt tau
required careful attention to siteSym text & some reorganization of modules
more work on special pos code
  • Property svn:eol-style set to native
  • Property svn:keywords set to Date Author Revision URL Id
File size: 156.1 KB
Line 
1# -*- coding: utf-8 -*-
2"""
3*GSASIIspc: Space group module*
4-------------------------------
5
6Space group interpretation routines. Note that space group information is
7stored in a :ref:`Space Group (SGData)<SGData_table>` object.
8
9"""
10########### SVN repository information ###################
11# $Date: 2015-08-07 16:54:12 +0000 (Fri, 07 Aug 2015) $
12# $Author: vondreele $
13# $Revision: 1951 $
14# $URL: trunk/GSASIIspc.py $
15# $Id: GSASIIspc.py 1951 2015-08-07 16:54:12Z vondreele $
16########### SVN repository information ###################
17import numpy as np
18import numpy.ma as ma
19import numpy.linalg as nl
20import scipy.optimize as so
21import math
22import sys
23import copy
24import os.path as ospath
25
26import GSASIIpath
27GSASIIpath.SetVersionNumber("$Revision: 1951 $")
28import pyspg
29
30npsind = lambda x: np.sin(x*np.pi/180.)
31npcosd = lambda x: np.cos(x*np.pi/180.)
32DEBUG = False
33   
34################################################################################
35#### Space group codes
36################################################################################
37
38def SpcGroup(SGSymbol):
39    """
40    Determines cell and symmetry information from a short H-M space group name
41
42    :param SGSymbol: space group symbol (string) with spaces between axial fields
43    :returns: (SGError,SGData)
44   
45       * SGError = 0 for no errors; >0 for errors (see SGErrors below for details)
46       * SGData - is a dict (see :ref:`Space Group object<SGData_table>`) with entries:
47       
48             * 'SpGrp': space group symbol, slightly cleaned up
49             * 'SGLaue':  one of '-1', '2/m', 'mmm', '4/m', '4/mmm', '3R',
50               '3mR', '3', '3m1', '31m', '6/m', '6/mmm', 'm3', 'm3m'
51             * 'SGInv': boolean; True if centrosymmetric, False if not
52             * 'SGLatt': one of 'P', 'A', 'B', 'C', 'I', 'F', 'R'
53             * 'SGUniq': one of 'a', 'b', 'c' if monoclinic, '' otherwise
54             * 'SGCen': cell centering vectors [0,0,0] at least
55             * 'SGOps': symmetry operations as [M,T] so that M*x+T = x'
56             * 'SGSys': one of 'triclinic', 'monoclinic', 'orthorhombic',
57               'tetragonal', 'rhombohedral', 'trigonal', 'hexagonal', 'cubic'
58             * 'SGPolax': one of ' ', 'x', 'y', 'x y', 'z', 'x z', 'y z',
59               'xyz', '111' for arbitrary axes
60             * 'SGPtGrp': one of 32 point group symbols (with some permutations), which
61                is filled by SGPtGroup, is external (KE) part of supersymmetry point group
62             * 'SSGKl': default internal (Kl) part of supersymmetry point group; modified
63                in supersymmetry stuff depending on chosen modulation vector for Mono & Ortho
64
65    """
66    LaueSym = ('-1','2/m','mmm','4/m','4/mmm','3R','3mR','3','3m1','31m','6/m','6/mmm','m3','m3m')
67    LattSym = ('P','A','B','C','I','F','R')
68    UniqSym = ('','','a','b','c','',)
69    SysSym = ('triclinic','monoclinic','orthorhombic','tetragonal','rhombohedral','trigonal','hexagonal','cubic')
70    SGData = {}
71    SGSymbol = SGSymbol.replace(':',' ')    #get rid of ':' in R space group symbols from some cif files
72    SGInfo = pyspg.sgforpy(SGSymbol)
73    SGData['SpGrp'] = SGSymbol.strip().lower().capitalize()
74    SGData['SGLaue'] = LaueSym[SGInfo[0]-1]
75    SGData['SGInv'] = bool(SGInfo[1])
76    SGData['SGLatt'] = LattSym[SGInfo[2]-1]
77    SGData['SGUniq'] = UniqSym[SGInfo[3]+1]
78    if SGData['SGLatt'] == 'P':
79        SGData['SGCen'] = np.array(([0,0,0],))
80    elif SGData['SGLatt'] == 'A':
81        SGData['SGCen'] = np.array(([0,0,0],[0,.5,.5]))
82    elif SGData['SGLatt'] == 'B':
83        SGData['SGCen'] = np.array(([0,0,0],[.5,0,.5]))
84    elif SGData['SGLatt'] == 'C':
85        SGData['SGCen'] = np.array(([0,0,0],[.5,.5,0,]))
86    elif SGData['SGLatt'] == 'I':
87        SGData['SGCen'] = np.array(([0,0,0],[.5,.5,.5]))
88    elif SGData['SGLatt'] == 'F':
89        SGData['SGCen'] = np.array(([0,0,0],[0,.5,.5],[.5,0,.5],[.5,.5,0,]))
90    elif SGData['SGLatt'] == 'R':
91        SGData['SGCen'] = np.array(([0,0,0],[1./3.,2./3.,2./3.],[2./3.,1./3.,1./3.]))
92    SGData['SGOps'] = []
93    for i in range(SGInfo[5]):
94        Mat = np.array(SGInfo[6][i])
95        Trns = np.array(SGInfo[7][i])
96        SGData['SGOps'].append([Mat,Trns])
97    if SGData['SGLaue'] in '-1':
98        SGData['SGSys'] = SysSym[0]
99    elif SGData['SGLaue'] in '2/m':
100        SGData['SGSys'] = SysSym[1]
101    elif SGData['SGLaue'] in 'mmm':
102        SGData['SGSys'] = SysSym[2]
103    elif SGData['SGLaue'] in ['4/m','4/mmm']:
104        SGData['SGSys'] = SysSym[3]
105    elif SGData['SGLaue'] in ['3R','3mR']:
106        SGData['SGSys'] = SysSym[4]
107    elif SGData['SGLaue'] in ['3','3m1','31m']:
108        SGData['SGSys'] = SysSym[5]
109    elif SGData['SGLaue'] in ['6/m','6/mmm']:
110        SGData['SGSys'] = SysSym[6]
111    elif SGData['SGLaue'] in ['m3','m3m']:
112        SGData['SGSys'] = SysSym[7]
113    SGData['SGPolax'] = SGpolar(SGData)
114    SGData['SGPtGrp'],SGData['SSGKl'] = SGPtGroup(SGData)
115    return SGInfo[8],SGData
116
117def SGErrors(IErr):
118    '''
119    Interprets the error message code from SpcGroup. Used in SpaceGroup.
120   
121    :param IErr: see SGError in :func:`SpcGroup`
122    :returns:
123        ErrString - a string with the error message or "Unknown error"
124    '''
125
126    ErrString = [' ',
127        'Less than 2 operator fields were found',
128        'Illegal Lattice type, not P, A, B, C, I, F or R',
129        'Rhombohedral lattice requires a 3-axis',
130        'Minus sign does not preceed 1, 2, 3, 4 or 6',
131        'Either a 5-axis anywhere or a 3-axis in field not allowed',
132        ' ',
133        'I for COMPUTED GO TO out of range.',
134        'An a-glide mirror normal to A not allowed',
135        'A b-glide mirror normal to B not allowed',
136        'A c-glide mirror normal to C not allowed',
137        'D-glide in a primitive lattice not allowed',
138        'A 4-axis not allowed in the 2nd operator field',
139        'A 6-axis not allowed in the 2nd operator field',
140        'More than 24 matrices needed to define group',
141        ' ',
142        'Improper construction of a rotation operator',
143        'Mirror following a / not allowed',
144        'A translation conflict between operators',
145        'The 2bar operator is not allowed',
146        '3 fields are legal only in R & m3 cubic groups',
147        'Syntax error. Expected I -4 3 d at this point',
148        ' ',
149        'A or B centered tetragonal not allowed',
150        ' ','unknown error in sgroup',' ',' ',' ',
151        'Illegal character in the space group symbol',
152        ]
153    try:
154        return ErrString[IErr]
155    except:
156        return "Unknown error"
157
158def SGpolar(SGData):
159    '''
160    Determine identity of polar axes if any
161    '''
162    POL = ('','x','y','x y','z','x z','y z','xyz','111')
163    NP = [1,2,4]
164    NPZ = [0,1]
165    for M,T in SGData['SGOps']:
166        for i in range(3):
167            if M[i][i] <= 0.: NP[i] = 0
168        if M[0][2] > 0: NPZ[0] = 8
169        if M[1][2] > 0: NPZ[1] = 0
170    NPol = (NP[0]+NP[1]+NP[2]+NPZ[0]*NPZ[1])*(1-int(SGData['SGInv']))
171    return POL[NPol]
172   
173def SGPtGroup(SGData):
174    '''
175    Determine point group of the space group - done after space group symbol has
176    been evaluated by SpcGroup. Only short symbols are allowed
177   
178    :param SGData: from :func SpcGroup
179    :returns: SSGPtGrp & SSGKl (only defaults for Mono & Ortho)
180    '''
181    Flds = SGData['SpGrp'].split()
182    if len(Flds) < 2:
183        return '',[]
184    if SGData['SGLaue'] == '-1':    #triclinic
185        if '-' in Flds[1]:
186            return '-1',[-1,]
187        else:
188            return '1',[1,]
189    elif SGData['SGLaue'] == '2/m': #monoclinic - default for 2D modulation vector
190        if '/' in SGData['SpGrp']:
191            return '2/m',[-1,1]
192        elif '2' in SGData['SpGrp']:
193            return '2',[-1,]
194        else:
195            return 'm',[1,]
196    elif SGData['SGLaue'] == 'mmm': #orthorhombic
197        if SGData['SpGrp'].count('2') == 3:
198            return '222',[-1,-1,-1]
199        elif SGData['SpGrp'].count('2') == 1:
200            if SGData['SGPolax'] == 'x':
201                return '2mm',[-1,1,1]
202            elif SGData['SGPolax'] == 'y':
203                return 'm2m',[1,-1,1]
204            elif SGData['SGPolax'] == 'z':
205                return 'mm2',[1,1,-1]
206        else:
207            return 'mmm',[1,1,1]
208    elif SGData['SGLaue'] == '4/m': #tetragonal
209        if '/' in SGData['SpGrp']:
210            return '4/m',[1,-1]
211        elif '-' in Flds[1]:
212            return '-4',[-1,]
213        else:
214            return '4',[1,]
215    elif SGData['SGLaue'] == '4/mmm':
216        if '/' in SGData['SpGrp']:
217            return '4/mmm',[1,-1,1,1]
218        elif '-' in Flds[1]:
219            if '2' in Flds[2]:
220                return '-42m',[-1,-1,1]
221            else:
222                return '-4m2',[-1,1,-1]             
223        elif '2' in Flds[2:]:
224            return '422',[1,-1,-1]
225        else:
226            return '4mm',[1,1,1]
227    elif SGData['SGLaue'] in ['3','3R']:  #trigonal/rhombohedral
228        if '-' in Flds[1]:
229            return '-3',[-1,]
230        else:
231            return '3',[1,]
232    elif SGData['SGLaue'] == '3mR' or 'R' in Flds[0]:
233        if '2' in Flds[2]:
234            return '32',[1,-1]
235        elif '-' in Flds[1]:
236            return '-3m',[-1,1]
237        else:
238            return '3m',[1,1]
239    elif SGData['SGLaue'] == '3m1':
240        if '2' in Flds[2]:
241            return '321',[1,-1,1]
242        elif '-' in Flds[1]:
243            return '-3m1',[-1,1,1]
244        else:
245            return '3m1',[1,1,1]
246    elif SGData['SGLaue'] == '31m':
247        if '2' in Flds[3]:
248            return '312',[1,1,-1]
249        elif '-' in Flds[1]:
250            return '-31m',[-1,1,1]
251        else:
252            return '31m',[1,1,1]
253    elif SGData['SGLaue'] == '6/m': #hexagonal
254        if '/' in SGData['SpGrp']:
255            return '6/m',[1,-1]
256        elif '-' in SGData['SpGrp']:
257            return '-6',[-1,]
258        else:
259            return '6',[1,]
260    elif SGData['SGLaue'] == '6/mmm':
261        if '/' in SGData['SpGrp']:
262            return '6/mmm',[1,-1,1,1]
263        elif '-' in Flds[1]:
264            if '2' in Flds[2]:
265                return '-62m',[-1,-1,1]
266            else:
267                return '-6m2',[-1,1,-1]                 
268        elif '2' in Flds[2:]:
269            return '622',[1,-1,-1]
270        else:
271            return '6mm',[1,1,1]   
272    elif SGData['SGLaue'] == 'm3':      #cubic - no (3+1) supersymmetry
273        if '2' in Flds[1]:
274            return '23',[]
275        else: 
276            return 'm3',[]
277    elif SGData['SGLaue'] == 'm3m':
278        if '4' in Flds[1]:
279            if '-' in Flds[1]:
280                return '-43m',[]
281            else:
282                return '432',[]
283        else:
284            return 'm-3m',[]
285   
286def SGPrint(SGData):
287    '''
288    Print the output of SpcGroup in a nicely formatted way. Used in SpaceGroup
289
290    :param SGData: from :func:`SpcGroup`
291    :returns:
292        SGText - list of strings with the space group details
293        SGTable - list of strings for each of the operations
294    '''
295    Mult = len(SGData['SGCen'])*len(SGData['SGOps'])*(int(SGData['SGInv'])+1)
296    SGText = []
297    SGText.append(' Space Group: '+SGData['SpGrp'])
298    CentStr = 'centrosymmetric'
299    if not SGData['SGInv']:
300        CentStr = 'non'+CentStr
301    if SGData['SGLatt'] in 'ABCIFR':
302        SGText.append(' The lattice is '+CentStr+' '+SGData['SGLatt']+'-centered '+SGData['SGSys'].lower())
303    else:
304        SGText.append(' The lattice is '+CentStr+' '+'primitive '+SGData['SGSys'].lower()) 
305    SGText.append(' The Laue symmetry is '+SGData['SGLaue'])
306    if 'SGPtGrp' in SGData:         #patch
307        SGText.append(' The lattice point group is '+SGData['SGPtGrp'])
308    SGText.append(' Multiplicity of a general site is '+str(Mult))
309    if SGData['SGUniq'] in ['a','b','c']:
310        SGText.append(' The unique monoclinic axis is '+SGData['SGUniq'])
311    if SGData['SGInv']:
312        SGText.append(' The inversion center is located at 0,0,0')
313    if SGData['SGPolax']:
314        SGText.append(' The location of the origin is arbitrary in '+SGData['SGPolax'])
315    SGText.append(' ')
316    if SGData['SGLatt'] == 'P':
317        SGText.append(' The equivalent positions are:\n')
318    else:   
319        SGText.append(' The equivalent positions are:')
320        SGText.append(' ('+Latt2text(SGData['SGLatt'])+')+\n')
321    SGTable = []
322    for i,Opr in enumerate(SGData['SGOps']):
323        SGTable.append('(%2d) %s'%(i+1,MT2text(Opr)))
324    return SGText,SGTable
325
326def AllOps(SGData):
327    '''
328    Returns a list of all operators for a space group, including those for
329    centering and a center of symmetry
330   
331    :param SGData: from :func:`SpcGroup`
332    :returns: (SGTextList,offsetList,symOpList,G2oprList) where
333
334      * SGTextList: a list of strings with formatted and normalized
335        symmetry operators.
336      * offsetList: a tuple of (dx,dy,dz) offsets that relate the GSAS-II
337        symmetry operation to the operator in SGTextList and symOpList.
338        these dx (etc.) values are added to the GSAS-II generated
339        positions to provide the positions that are generated
340        by the normalized symmetry operators.       
341      * symOpList: a list of tuples with the normalized symmetry
342        operations as (M,T) values
343        (see ``SGOps`` in the :ref:`Space Group object<SGData_table>`)
344      * G2oprList: The GSAS-II operations for each symmetry operation as
345        a tuple with (center,mult,opnum), where center is (0,0,0), (0.5,0,0),
346        (0.5,0.5,0.5),...; where mult is 1 or -1 for the center of symmetry
347        and opnum is the number for the symmetry operation, in ``SGOps``
348        (starting with 0).
349    '''
350    SGTextList = []
351    offsetList = []
352    symOpList = []
353    G2oprList = []
354    onebar = (1,)
355    if SGData['SGInv']:
356        onebar += (-1,)
357    for cen in SGData['SGCen']:
358        for mult in onebar:
359            for j,(M,T) in enumerate(SGData['SGOps']):
360                offset = [0,0,0]
361                Tprime = (mult*T)+cen
362                for i in range(3):
363                    while Tprime[i] < 0:
364                        Tprime[i] += 1
365                        offset[i] += 1
366                    while Tprime[i] >= 1:
367                        Tprime[i] += -1
368                        offset[i] += -1
369                Opr = [mult*M,Tprime]
370                OPtxt = MT2text(Opr)
371                SGTextList.append(OPtxt.replace(' ',''))
372                offsetList.append(tuple(offset))
373                symOpList.append((mult*M,Tprime))
374                G2oprList.append((cen,mult,j))
375    return SGTextList,offsetList,symOpList,G2oprList
376   
377def MT2text(Opr):
378    "From space group matrix/translation operator returns text version"
379    XYZ = ('-Z','-Y','-X','X-Y','ERR','Y-X','X','Y','Z')
380    TRA = ('   ','ERR','1/6','1/4','1/3','ERR','1/2','ERR','2/3','3/4','5/6','ERR')
381    Fld = ''
382    M,T = Opr
383    for j in range(3):
384        IJ = int(round(2*M[j][0]+3*M[j][1]+4*M[j][2]+4))%12
385        IK = int(round(T[j]*12))%12
386        if IK:
387            if IJ < 3:
388                Fld += (TRA[IK]+XYZ[IJ]).rjust(5)
389            else:
390                Fld += (TRA[IK]+'+'+XYZ[IJ]).rjust(5)
391        else:
392            Fld += XYZ[IJ].rjust(5)
393        if j != 2: Fld += ', '
394    return Fld
395   
396def Latt2text(Latt):
397    "From lattice type ('P',A', etc.) returns ';' delimited cell centering vectors"
398    lattTxt = {'A':'0,0,0; 0,1/2,1/2','B':'0,0,0; 1/2,0,1/2',
399        'C':'0,0,0; 1/2,1/2,0','I':'0,0,0; 1/2,1/2,1/2',
400        'F':'0,0,0; 0,1/2,1/2; 1/2,0,1/2; 1/2,1/2,0',
401        'R':'0,0,0; 1/3,2/3,2/3; 2/3,1/3,1/3','P':'0,0,0'}
402    return lattTxt[Latt]   
403       
404def SpaceGroup(SGSymbol):
405    '''
406    Print the output of SpcGroup in a nicely formatted way.
407
408    :param SGSymbol: space group symbol (string) with spaces between axial fields
409    :returns: nothing
410    '''
411    E,A = SpcGroup(SGSymbol)
412    if E > 0:
413        print SGErrors(E)
414        return
415    for l in SGPrint(A):
416        print l
417       
418################################################################################
419#### Superspace group codes
420################################################################################
421       
422def SSpcGroup(SGData,SSymbol):
423    """
424    Determines supersymmetry information from superspace group name; currently only for (3+1) superlattices
425
426    :param SGData: space group data structure as defined in SpcGroup above (see :ref:`SGData<SGData_table>`).
427    :param SSymbol: superspace group symbol extension (string) defining modulation direction & generator info.
428    :returns: (SSGError,SSGData)
429   
430       * SGError = 0 for no errors; >0 for errors (see SGErrors below for details)
431       * SSGData - is a dict (see :ref:`Superspace Group object<SSGData_table>`) with entries:
432       
433             * 'SSpGrp': superspace group symbol extension to space group symbol, accidental spaces removed
434             * 'SSGCen': 4D cell centering vectors [0,0,0,0] at least
435             * 'SSGOps': 4D symmetry operations as [M,T] so that M*x+T = x'
436
437    """
438   
439    def checkModSym():
440        '''
441        Checks to see if proposed modulation form is allowed for Laue group
442        '''
443        if LaueId in [0,] and LaueModId in [0,]:
444            return True
445        elif LaueId in [1,]:
446            try:
447                if modsym.index('1/2') != ['A','B','C'].index(SGData['SGLatt']):
448                    return False
449                if 'I'.index(SGData['SGLatt']) and modsym.count('1/2') not in [0,2]:
450                    return False
451            except ValueError:
452                pass
453            if SGData['SGUniq'] == 'a' and LaueModId in [5,6,7,8,9,10,]:
454                return True
455            elif SGData['SGUniq'] == 'b' and LaueModId in [3,4,13,14,15,16,]:
456                return True
457            elif SGData['SGUniq'] == 'c' and LaueModId in [1,2,19,20,21,22,]:
458                return True
459        elif LaueId in [2,] and LaueModId in [i+7 for i in range(18)]:
460            try:
461                if modsym.index('1/2') != ['A','B','C'].index(SGData['SGLatt']):
462                    return False
463                if SGData['SGLatt'] in ['I','F',] and modsym.index('1/2'):
464                    return False
465            except ValueError:
466                pass
467            return True
468        elif LaueId in [3,4,] and LaueModId in [19,22,]:
469            try:
470                if SGData['SGLatt'] == 'I' and modsym.count('1/2'):
471                    return False
472            except ValueError:
473                pass
474            return True
475        elif LaueId in [7,8,9,] and LaueModId in [19,25,]:
476            if (SGData['SGLatt'] == 'R' or SGData['SGPtGrp'] in ['3m1','-3m1']) and modsym.count('1/3'):
477                return False
478            return True
479        elif LaueId in [10,11,] and LaueModId in [19,]:
480            return True
481        return False
482       
483    def fixMonoOrtho():
484        mod = ''.join(modsym).replace('1/2','0').replace('1','0')
485        if SGData['SGPtGrp'] in ['2','m']:  #OK
486            if mod in ['a00','0b0','00g']:
487                result = [i*-1 for i in SGData['SSGKl']]
488            else:
489                result = SGData['SSGKl'][:]
490            if '/' in mod:
491                return [i*-1 for i in result]
492            else:
493                return result
494        elif SGData['SGPtGrp'] == '2/m':    #OK
495            if mod in ['a00','0b0','00g']:
496                result =  SGData['SSGKl'][:]
497            else:
498                result = [i*-1 for i in SGData['SSGKl']]
499            if '/' in mod:
500                return [i*-1 for i in result]
501            else:
502                return result
503        else:   #orthorhombic
504            return [-SSGKl[i] if mod[i] in ['a','b','g'] else SSGKl[i] for i in range(3)]
505               
506    def extendSSGOps(SSGOps):
507        nOps = len(SSGOps)
508        for OpA in SSGOps:
509            OpAtxt = SSMT2text(OpA)
510            if 't' not in OpAtxt:
511                continue
512            for OpB in SSGOps:
513                OpBtxt = SSMT2text(OpB)
514                if 't' not in OpBtxt:
515                    continue
516                OpC = list(SGProd(OpB,OpA))
517                OpC[1] %= 1.
518                OpCtxt = SSMT2text(OpC)
519#                print OpAtxt.replace(' ','')+' * '+OpBtxt.replace(' ','')+' = '+OpCtxt.replace(' ','')
520                for k,OpD in enumerate(SSGOps):
521                    OpDtxt = SSMT2text(OpD)
522                    if 't' in OpDtxt:
523                        continue
524#                    print '    ('+OpCtxt.replace(' ','')+' = ? '+OpDtxt.replace(' ','')+')'
525                    if OpCtxt == OpDtxt:
526                        continue
527                    elif OpCtxt.split(',')[:3] == OpDtxt.split(',')[:3]:
528                        if 't' not in OpDtxt:
529                            SSGOps[k] = OpC
530#                            print k,'   new:',OpCtxt.replace(' ','')
531                            break
532                        else:
533                            OpCtxt = OpCtxt.replace(' ','')
534                            OpDtxt = OpDtxt.replace(' ','')
535                            Txt = OpCtxt+' conflict with '+OpDtxt
536                            print Txt
537                            return False,Txt
538        return True,SSGOps
539       
540    def findMod(modSym):
541        for a in ['a','b','g']:
542            if a in modSym:
543                return a
544               
545    def genSSGOps():
546        SSGOps = SSGData['SSGOps'][:]
547        iFrac = {}
548        for i,frac in enumerate(SSGData['modSymb']):
549            if frac in ['1/2','1/3','1/4','1/6','1']:
550                iFrac[i] = frac+'.'
551#        print SGData['SpGrp']+SSymbol
552#        print 'SSGKl',SSGKl,'genQ',genQ,'iFrac',iFrac,'modSymb',SSGData['modSymb']
553# set identity & 1,-1; triclinic
554        SSGOps[0][0][3,3] = 1.
555## expand if centrosymmetric
556#        if SGData['SGInv']:
557#            SSGOps += [[-1*M,V] for M,V in SSGOps[:]]
558# monoclinic - all done & all checked
559        if SGData['SGPtGrp'] in ['2','m']:  #OK
560            SSGOps[1][0][3,3] = SSGKl[0]
561            SSGOps[1][1][3] = genQ[0]
562            for i in iFrac:
563                SSGOps[1][0][3,i] = -SSGKl[0]
564        elif SGData['SGPtGrp'] == '2/m':    #OK
565            SSGOps[1][0][3,3] = SSGKl[1]
566            if gensym:
567                SSGOps[1][1][3] = 0.5
568            for i in iFrac:
569                SSGOps[1][0][3,i] = SSGKl[0]
570           
571# orthorhombic - all OK not fully checked
572        elif SGData['SGPtGrp'] in ['222','mm2','m2m','2mm']:    #OK
573            if SGData['SGPtGrp'] == '222':
574                OrOps = {'g':{0:[1,3],1:[2,3]},'a':{1:[1,2],2:[1,3]},'b':{2:[3,2],0:[1,2]}} #OK
575            elif SGData['SGPtGrp'] == 'mm2':
576                OrOps = {'g':{0:[1,3],1:[2,3]},'a':{1:[2,1],2:[3,1]},'b':{0:[1,2],2:[3,2]}} #OK
577            elif SGData['SGPtGrp'] == 'm2m':
578                OrOps = {'b':{0:[1,2],2:[3,2]},'g':{0:[1,3],1:[2,3]},'a':{1:[2,1],2:[3,1]}} #OK
579            elif SGData['SGPtGrp'] == '2mm':
580                OrOps = {'a':{1:[2,1],2:[3,1]},'b':{0:[1,2],2:[3,2]},'g':{0:[1,3],1:[2,3]}} #OK
581            a = findMod(SSGData['modSymb'])
582            OrFrac = OrOps[a]
583            for j in iFrac:
584                for i in OrFrac[j]:
585                    SSGOps[i][0][3,j] = -2.*eval(iFrac[j])*SSGKl[i-1]
586            for i in [0,1,2]:
587                SSGOps[i+1][0][3,3] = SSGKl[i]
588                SSGOps[i+1][1][3] = genQ[i]
589                E,SSGOps = extendSSGOps(SSGOps)
590                if not E:
591                    return E,SSGOps
592        elif SGData['SGPtGrp'] == 'mmm':    #OK
593            OrOps = {'g':{0:[1,3],1:[2,3]},'a':{1:[2,1],2:[3,1]},'b':{0:[1,2],2:[3,2]}} 
594            a = findMod(SSGData['modSymb'])
595            if a == 'g':
596                SSkl = [1,1,1]
597            elif a == 'a':
598                SSkl = [-1,1,-1]
599            else:
600                SSkl = [1,-1,-1]
601            OrFrac = OrOps[a]
602            for j in iFrac:
603                for i in OrFrac[j]:
604                    SSGOps[i][0][3,j] = -2.*eval(iFrac[j])*SSkl[i-1]
605            for i in [0,1,2]:
606                SSGOps[i+1][0][3,3] = SSkl[i]
607                SSGOps[i+1][1][3] = genQ[i]
608                E,SSGOps = extendSSGOps(SSGOps)
609                if not E:
610                    return E,SSGOps               
611# tetragonal - all done & checked
612        elif SGData['SGPtGrp'] == '4':  #OK
613            SSGOps[1][0][3,3] = SSGKl[0]
614            SSGOps[1][1][3] = genQ[0]
615            if '1/2' in SSGData['modSymb']:
616                SSGOps[1][0][3,1] = -1
617        elif SGData['SGPtGrp'] == '-4': #OK
618            SSGOps[1][0][3,3] = SSGKl[0]
619            if '1/2' in SSGData['modSymb']:
620                SSGOps[1][0][3,1] = 1
621        elif SGData['SGPtGrp'] in ['4/m',]: #OK
622            if '1/2' in SSGData['modSymb']:
623                SSGOps[1][0][3,1] = -SSGKl[0]
624            for i,j in enumerate([1,3]):
625                SSGOps[j][0][3,3] = 1
626                if genQ[i]:
627                    SSGOps[j][1][3] = genQ[i]
628                E,SSGOps = extendSSGOps(SSGOps)
629                if not E:
630                    return E,SSGOps
631        elif SGData['SGPtGrp'] in ['422','4mm','-42m','-4m2',]: #OK
632            iGens = [1,4,5]
633            if SGData['SGPtGrp'] in ['4mm','-4m2',]:
634                iGens = [1,6,7]
635            for i,j in enumerate(iGens):
636                if '1/2' in SSGData['modSymb'] and i < 2:
637                    SSGOps[j][0][3,1] = SSGKl[i]
638                SSGOps[j][0][3,3] = SSGKl[i]
639                if genQ[i]:
640                    if 's' in gensym and j == 6:
641                        SSGOps[j][1][3] = -genQ[i]
642                    else:
643                        SSGOps[j][1][3] = genQ[i]
644                E,SSGOps = extendSSGOps(SSGOps)
645                if not E:
646                    return E,SSGOps
647        elif SGData['SGPtGrp'] in ['4/mmm',]:#OK
648            if '1/2' in SSGData['modSymb']:
649                SSGOps[1][0][3,1] = -SSGKl[0]
650                SSGOps[6][0][3,1] = SSGKl[1]
651                if modsym:
652                   SSGOps[1][1][3]  = -genQ[3]
653            for i,j in enumerate([1,2,6,7]):
654                SSGOps[j][0][3,3] = 1
655                SSGOps[j][1][3] = genQ[i]
656                E,Result = extendSSGOps(SSGOps)
657                if not E:
658                    return E,Result
659                else:
660                    SSGOps = Result
661               
662# trigonal - all done & checked
663        elif SGData['SGPtGrp'] == '3':  #OK
664            SSGOps[1][0][3,3] = SSGKl[0]
665            if '1/3' in SSGData['modSymb']:
666                SSGOps[1][0][3,1] = -1
667            SSGOps[1][1][3] = genQ[0]
668        elif SGData['SGPtGrp'] == '-3': #OK
669            SSGOps[1][0][3,3] = -SSGKl[0]
670            if '1/3' in SSGData['modSymb']:
671                SSGOps[1][0][3,1] = -1
672            SSGOps[1][1][3] = genQ[0]
673        elif SGData['SGPtGrp'] in ['312','3m','-3m','-3m1','3m1']:   #OK
674            if '1/3' in SSGData['modSymb']:
675                SSGOps[1][0][3,1] = -1
676            for i,j in enumerate([1,5]):
677                if SGData['SGPtGrp'] in ['3m','-3m']:
678                    SSGOps[j][0][3,3] = 1
679                else:                   
680                    SSGOps[j][0][3,3] = SSGKl[i+1]
681                if genQ[i]:
682                    SSGOps[j][1][3] = genQ[i]
683        elif SGData['SGPtGrp'] in ['321','32']:   #OK
684            for i,j in enumerate([1,4]):
685                SSGOps[j][0][3,3] = SSGKl[i]
686                if genQ[i]:
687                    SSGOps[j][1][3] = genQ[i]
688        elif SGData['SGPtGrp'] in ['31m','-31m']:   #OK
689            ids = [1,3]
690            if SGData['SGPtGrp'] == '-31m':
691                ids = [1,3]
692            if '1/3' in SSGData['modSymb']:
693                SSGOps[ids[0]][0][3,1] = -SSGKl[0]
694            for i,j in enumerate(ids):
695                SSGOps[j][0][3,3] = 1
696                if genQ[i+1]:
697                    SSGOps[j][1][3] = genQ[i+1]
698                     
699# hexagonal all done & checked
700        elif SGData['SGPtGrp'] == '6':  #OK
701            SSGOps[1][0][3,3] = SSGKl[0]
702            SSGOps[1][1][3] = genQ[0]
703        elif SGData['SGPtGrp'] == '-6': #OK
704            SSGOps[1][0][3,3] = SSGKl[0]
705        elif SGData['SGPtGrp'] in ['6/m',]: #OK
706            SSGOps[1][0][3,3] = -SSGKl[1]
707            SSGOps[1][1][3] = genQ[0]
708            SSGOps[2][1][3] = genQ[1]
709        elif SGData['SGPtGrp'] in ['622',]: #OK
710            for i,j in enumerate([1,8,9]):
711                SSGOps[j][0][3,3] = SSGKl[i]
712                if genQ[i]:
713                    SSGOps[j][1][3] = genQ[i]
714                E,SSGOps = extendSSGOps(SSGOps)
715           
716        elif SGData['SGPtGrp'] in ['6mm','-62m','-6m2',]: #OK
717            for i,j in enumerate([1,6,7]):
718                SSGOps[j][0][3,3] = SSGKl[i]
719                if genQ[i]:
720                    SSGOps[j][1][3] = genQ[i]
721                E,SSGOps = extendSSGOps(SSGOps)
722        elif SGData['SGPtGrp'] in ['6/mmm',]: # OK
723            for i,j in enumerate([1,2,10,11]):
724                SSGOps[j][0][3,3] = 1
725                if genQ[i]:
726                    SSGOps[j][1][3] = genQ[i]
727                E,SSGOps = extendSSGOps(SSGOps)
728        elif SGData['SGPtGrp'] in ['1','-1']: #triclinic - done
729            return True,SSGOps
730        E,SSGOps = extendSSGOps(SSGOps)
731        return E,SSGOps
732       
733    def specialGen(gensym,modsym):
734        sym = ''.join(gensym)
735        if SGData['SGPtGrp'] in ['2/m',] and 'n' in SGData['SpGrp']:
736            if 's' in sym:
737                gensym = 'ss'
738        if SGData['SGPtGrp'] in ['-62m',] and sym == '00s':
739            gensym = '0ss'
740        elif SGData['SGPtGrp'] in ['222',]:
741            if sym == '00s':
742                gensym = '0ss'
743            elif sym == '0s0':
744                gensym = 'ss0'
745            elif sym == 's00':
746                gensym = 's0s'
747        elif SGData['SGPtGrp'] in ['mmm',]:
748            if 'g' in modsym:
749                if sym == 's00':
750                    gensym = 's0s'
751                elif sym == '0s0':
752                    gensym = '0ss'
753            elif 'a' in modsym:
754                if sym == '0s0':
755                    gensym = 'ss0'
756                elif sym == '00s':
757                    gensym = 's0s'
758            elif 'b' in modsym:
759                if sym == '00s':
760                    gensym = '0ss'
761                elif sym == 's00':
762                    gensym = 'ss0'
763        return gensym
764                   
765    def checkGen(gensym):
766        '''
767    GenSymList = ['','s','0s','s0', '00s','0s0','s00','s0s','ss0','0ss','q00','0q0','00q','qq0','q0q', '0qq',
768        'q','qqs','s0s0','00ss','s00s','t','t00','t0','h','h00','000s']
769        '''
770        sym = ''.join(gensym)
771# monoclinic - all done
772        if str(SSGKl) == '[-1]' and sym == 's':
773            return False
774        elif SGData['SGPtGrp'] in ['2/m',]:
775            if str(SSGKl) == '[-1, 1]' and sym == '0s':
776                return False
777            elif str(SSGKl) == '[1, -1]' and sym == 's0':
778                return False
779#orthorhombic - all
780        elif SGData['SGPtGrp'] in ['222',] and sym not in ['','s00','0s0','00s']:
781            return False 
782        elif SGData['SGPtGrp'] in ['2mm','m2m','mm2','mmm'] and sym not in ['',]+GenSymList[4:16]:
783            return False 
784#tetragonal - all done
785        elif SGData['SGPtGrp'] in ['4',] and sym not in ['','s','q']:
786            return False 
787        elif SGData['SGPtGrp'] in ['-4',] and sym not in ['',]:
788            return False             
789        elif SGData['SGPtGrp'] in ['4/m',] and sym not in ['','s0','q0']:
790            return False
791        elif SGData['SGPtGrp'] in ['422',] and sym not in ['','q00','s00']:
792            return False         
793        elif SGData['SGPtGrp'] in ['4mm',] and sym not in ['','ss0','s0s','0ss','00s','qq0','qqs']:
794            return False
795        elif SGData['SGPtGrp'] in ['-4m2',] and sym not in ['','0s0','0q0']:
796            return False
797        elif SGData['SGPtGrp'] in ['-42m',] and sym not in ['','0ss','00q',]:
798            return False
799        elif SGData['SGPtGrp'] in ['4/mmm',] and sym not in ['','s00s','s0s0','00ss','000s',]:
800            return False
801#trigonal/rhombohedral - all done
802        elif SGData['SGPtGrp'] in ['3',] and sym not in ['','t']:
803            return False 
804        elif SGData['SGPtGrp'] in ['-3',] and sym not in ['',]:
805            return False 
806        elif SGData['SGPtGrp'] in ['32',] and sym not in ['','t0']:
807            return False 
808        elif SGData['SGPtGrp'] in ['321','312'] and sym not in ['','t00']:
809            return False 
810        elif SGData['SGPtGrp'] in ['3m','-3m'] and sym not in ['','0s']:
811            return False 
812        elif SGData['SGPtGrp'] in ['3m1','-3m1'] and sym not in ['','0s0']:
813            return False 
814        elif SGData['SGPtGrp'] in ['31m','-31m'] and sym not in ['','00s']:
815            return False 
816#hexagonal - all done
817        elif SGData['SGPtGrp'] in ['6',] and sym not in ['','s','h','t']:
818            return False 
819        elif SGData['SGPtGrp'] in ['-6',] and sym not in ['',]:
820            return False
821        elif SGData['SGPtGrp'] in ['6/m',] and sym not in ['','s0']:
822            return False
823        elif SGData['SGPtGrp'] in ['622',] and sym not in ['','h00','t00','s00']:
824            return False         
825        elif SGData['SGPtGrp'] in ['6mm',] and sym not in ['','ss0','s0s','0ss']:
826            return False
827        elif SGData['SGPtGrp'] in ['-6m2',] and sym not in ['','0s0']:
828            return False
829        elif SGData['SGPtGrp'] in ['-62m',] and sym not in ['','00s']:
830            return False
831        elif SGData['SGPtGrp'] in ['6/mmm',] and sym not in ['','s00s','s0s0','00ss']:
832            return False
833        return True
834       
835    LaueModList = [
836        'abg','ab0','ab1/2','a0g','a1/2g',  '0bg','1/2bg','a00','a01/2','a1/20',
837        'a1/21/2','a01','a10','0b0','0b1/2', '1/2b0','1/2b1/2','0b1','1b0','00g',
838        '01/2g','1/20g','1/21/2g','01g','10g', '1/31/3g']
839    LaueList = ['-1','2/m','mmm','4/m','4/mmm','3R','3mR','3','3m1','31m','6/m','6/mmm','m3','m3m']
840    GenSymList = ['','s','0s','s0', '00s','0s0','s00','s0s','ss0','0ss','q00','0q0','00q','qq0','q0q', '0qq',
841        'q','qqs','s0s0','00ss','s00s','t','t00','t0','h','h00','000s']
842    Fracs = {'1/2':0.5,'1/3':1./3,'1':1.0,'0':0.,'s':.5,'t':1./3,'q':.25,'h':1./6,'a':0.,'b':0.,'g':0.}
843    LaueId = LaueList.index(SGData['SGLaue'])
844    if SGData['SGLaue'] in ['m3','m3m']:
845        return '(3+1) superlattices not defined for cubic space groups',None
846    elif SGData['SGLaue'] in ['3R','3mR']:
847        return '(3+1) superlattices not defined for rhombohedral settings - use hexagonal setting',None
848    try:
849        modsym,gensym = splitSSsym(SSymbol)
850    except ValueError:
851        return 'Error in superspace symbol '+SSymbol,None
852    if ''.join(gensym) not in GenSymList:
853        return 'unknown generator symbol '+''.join(gensym),None
854    try:
855        LaueModId = LaueModList.index(''.join(modsym))
856    except ValueError:
857        return 'Unknown modulation symbol '+''.join(modsym),None
858    if not checkModSym():
859        return 'Modulation '+''.join(modsym)+' not consistent with space group '+SGData['SpGrp'],None
860    modQ = [Fracs[mod] for mod in modsym]
861    SSGKl = SGData['SSGKl'][:]
862    if SGData['SGLaue'] in ['2/m','mmm']:
863        SSGKl = fixMonoOrtho()
864    if len(gensym) and len(gensym) != len(SSGKl):
865        return 'Wrong number of items in generator symbol '+''.join(gensym),None
866    if not checkGen(gensym):
867        return 'Generator '+''.join(gensym)+' not consistent with space group '+SGData['SpGrp'],None
868    gensym = specialGen(gensym,modsym)
869    genQ = [Fracs[mod] for mod in gensym]
870    if not genQ:
871        genQ = [0,0,0,0]
872    SSGData = {'SSpGrp':SGData['SpGrp']+SSymbol,'modQ':modQ,'modSymb':modsym,'SSGKl':SSGKl}
873    SSCen = np.zeros((len(SGData['SGCen']),4))
874    for icen,cen in enumerate(SGData['SGCen']):
875        SSCen[icen,0:3] = cen
876    SSCen[0] = np.zeros(4)
877    SSGData['SSGCen'] = SSCen
878    SSGData['SSGOps'] = []
879    for iop,op in enumerate(SGData['SGOps']):
880        T = np.zeros(4)
881        ssop = np.zeros((4,4))
882        ssop[:3,:3] = op[0]
883        T[:3] = op[1]
884        SSGData['SSGOps'].append([ssop,T])
885    E,Result = genSSGOps()
886    if E:
887        SSGData['SSGOps'] = Result
888        if DEBUG:
889            print 'Super spacegroup operators for '+SSGData['SSpGrp']
890            for Op in Result:
891                print SSMT2text(Op).replace(' ','')
892            if SGData['SGInv']:                                 
893                for Op in Result:
894                    Op = [-Op[0],-Op[1]%1.]
895                    print SSMT2text(Op).replace(' ','')                                 
896        return None,SSGData
897    else:
898        return Result+'\nOperator conflict - incorrect superspace symbol',None
899
900def splitSSsym(SSymbol):
901    '''
902    Splits supersymmetry symbol into two lists of strings
903    '''
904    modsym,gensym = SSymbol.replace(' ','').split(')')
905    nfrac = modsym.count('/')
906    modsym = modsym.lstrip('(')
907    if nfrac == 0:
908        modsym = list(modsym)
909    elif nfrac == 1:
910        pos = modsym.find('/')
911        if pos == 1:
912            modsym = [modsym[:3],modsym[3],modsym[4]]
913        elif pos == 2:
914            modsym = [modsym[0],modsym[1:4],modsym[4]]
915        else:
916            modsym = [modsym[0],modsym[1],modsym[2:]]
917    else:
918        lpos = modsym.find('/')
919        rpos = modsym.rfind('/')
920        if lpos == 1 and rpos == 4:
921            modsym = [modsym[:3],modsym[3:6],modsym[6]]
922        elif lpos == 1 and rpos == 5:
923            modsym = [modsym[:3],modsym[3],modsym[4:]]
924        else:
925            modsym = [modsym[0],modsym[1:4],modsym[4:]]
926    gensym = list(gensym)
927    return modsym,gensym
928       
929def SSGPrint(SGData,SSGData):
930    '''
931    Print the output of SSpcGroup in a nicely formatted way. Used in SSpaceGroup
932
933    :param SGData: space group data structure as defined in SpcGroup above.
934    :param SSGData: from :func:`SSpcGroup`
935    :returns:
936        SSGText - list of strings with the superspace group details
937        SGTable - list of strings for each of the operations
938    '''
939    Mult = len(SSGData['SSGCen'])*len(SSGData['SSGOps'])*(int(SGData['SGInv'])+1)
940    SSGText = []
941    SSGText.append(' Superspace Group: '+SSGData['SSpGrp'])
942    CentStr = 'centrosymmetric'
943    if not SGData['SGInv']:
944        CentStr = 'non'+CentStr
945    if SGData['SGLatt'] in 'ABCIFR':
946        SSGText.append(' The lattice is '+CentStr+' '+SGData['SGLatt']+'-centered '+SGData['SGSys'].lower())
947    else:
948        SSGText.append(' The superlattice is '+CentStr+' '+'primitive '+SGData['SGSys'].lower())       
949    SSGText.append(' The Laue symmetry is '+SGData['SGLaue'])
950    SSGText.append(' The superlattice point group is '+SGData['SGPtGrp']+', '+''.join([str(i) for i in SSGData['SSGKl']]))
951    SSGText.append(' The number of superspace group generators is '+str(len(SGData['SSGKl'])))
952    SSGText.append(' Multiplicity of a general site is '+str(Mult))
953    if SGData['SGUniq'] in ['a','b','c']:
954        SSGText.append(' The unique monoclinic axis is '+SGData['SGUniq'])
955    if SGData['SGInv']:
956        SSGText.append(' The inversion center is located at 0,0,0')
957    if SGData['SGPolax']:
958        SSGText.append(' The location of the origin is arbitrary in '+SGData['SGPolax'])
959    SSGText.append(' ')
960    if len(SSGData['SSGCen']) > 1:
961        SSGText.append(' The equivalent positions are:')
962        SSGText.append(' ('+SSLatt2text(SSGData['SSGCen'])+')+\n')
963    else:
964        SSGText.append(' The equivalent positions are:\n')
965    SSGTable = []
966    for i,Opr in enumerate(SSGData['SSGOps']):
967        SSGTable.append('(%2d) %s'%(i+1,SSMT2text(Opr)))
968    return SSGText,SSGTable
969   
970def SSGModCheck(Vec,modSymb):
971    ''' Checks modulation vector compatibility with supersymmetry space group symbol.
972    Superspace group symbol takes precidence & the vector will be modified accordingly
973    '''
974    Fracs = {'1/2':0.5,'1/3':1./3,'1':1.0,'0':0.,'a':0.,'b':0.,'g':0.}
975    modQ = [Fracs[mod] for mod in modSymb]
976    Vec = [0.1 if (vec == 0.0 and mod in ['a','b','g']) else vec for [vec,mod] in zip(Vec,modSymb)]
977    return [Q if mod not in ['a','b','g'] and vec != Q else vec for [vec,mod,Q] in zip(Vec,modSymb,modQ)],  \
978        [True if mod in ['a','b','g'] else False for mod in modSymb]
979
980def SSMT2text(Opr):
981    "From superspace group matrix/translation operator returns text version"
982    XYZS = ('x','y','z','t')    #Stokes, Campbell & van Smaalen notation
983    TRA = ('   ','ERR','1/6','1/4','1/3','ERR','1/2','ERR','2/3','3/4','5/6','ERR')
984    Fld = ''
985    M,T = Opr
986    for j in range(4):
987        IJ = ''
988        for k in range(4):
989            txt = str(int(round(M[j][k])))
990            txt = txt.replace('1',XYZS[k]).replace('0','')
991            if '2' in txt:
992                txt += XYZS[k]
993            if IJ and M[j][k] > 0:
994                IJ += '+'+txt
995            else:
996                IJ += txt
997        IK = int(round(T[j]*12))%12
998        if IK:
999            if not IJ:
1000                break
1001            if IJ[0] == '-':
1002                Fld += (TRA[IK]+IJ).rjust(8)
1003            else:
1004                Fld += (TRA[IK]+'+'+IJ).rjust(8)
1005        else:
1006            Fld += IJ.rjust(8)
1007        if j != 3: Fld += ', '
1008    return Fld
1009   
1010def SSLatt2text(SSGCen):
1011    "Lattice centering vectors to text"
1012    lattTxt = ''
1013    for vec in SSGCen:
1014        lattTxt += ' '
1015        for item in vec:
1016            if int(item*12.):
1017                lattTxt += '1/%d,'%(12/int(item*12))
1018            else:
1019                lattTxt += '0,'
1020        lattTxt = lattTxt.rstrip(',')
1021        lattTxt += ';'
1022    lattTxt = lattTxt.rstrip(';').lstrip(' ')
1023    return lattTxt
1024       
1025def SSpaceGroup(SGSymbol,SSymbol):
1026    '''
1027    Print the output of SSpcGroup in a nicely formatted way.
1028
1029    :param SGSymbol: space group symbol with spaces between axial fields.
1030    :param SSymbol: superspace group symbol extension (string).
1031    :returns: nothing
1032    '''
1033
1034    E,A = SpcGroup(SGSymbol)
1035    if E > 0:
1036        print SGErrors(E)
1037        return
1038    E,B = SSpcGroup(A,SSymbol)   
1039    if E > 0:
1040        print E
1041        return
1042    for l in SSGPrint(B):
1043        print l
1044       
1045def SGProd(OpA,OpB):
1046    '''
1047    Form space group operator product. OpA & OpB are [M,V] pairs;
1048        both must be of same dimension (3 or 4). Returns [M,V] pair
1049    '''
1050    A,U = OpA
1051    B,V = OpB
1052    M = np.inner(B,A.T)
1053    W = np.inner(B,U)+V
1054    return M,W
1055       
1056def MoveToUnitCell(xyz):
1057    '''
1058    Translates a set of coordinates so that all values are >=0 and < 1
1059
1060    :param xyz: a list or numpy array of fractional coordinates
1061    :returns: XYZ - numpy array of new coordinates now 0 or greater and less than 1
1062    '''
1063    XYZ = np.array([(x-int(x))%1.0 for x in xyz])
1064    cell = np.asarray(np.rint(xyz-XYZ),dtype=np.int32)
1065    return XYZ,cell
1066       
1067def Opposite(XYZ,toler=0.0002):
1068    '''
1069    Gives opposite corner, edge or face of unit cell for position within tolerance.
1070        Result may be just outside the cell within tolerance
1071
1072    :param XYZ: 0 >= np.array[x,y,z] > 1 as by MoveToUnitCell
1073    :param toler: unit cell fraction tolerance making opposite
1074    :returns:
1075        XYZ: dict of opposite positions; key=unit cell & always contains XYZ
1076    '''
1077    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]]
1078    TB = np.where(abs(XYZ-1)<toler,-1,0)+np.where(abs(XYZ)<toler,1,0)
1079    perm = TB*perm3
1080    cperm = ['%d,%d,%d'%(i,j,k) for i,j,k in perm]
1081    D = dict(zip(cperm,perm))
1082    new = {}
1083    for key in D:
1084        new[key] = np.array(D[key])+np.array(XYZ)
1085    return new
1086       
1087def GenAtom(XYZ,SGData,All=False,Uij=[],Move=True):
1088    '''
1089    Generates the equivalent positions for a specified coordinate and space group
1090
1091    :param XYZ: an array, tuple or list containing 3 elements: x, y & z
1092    :param SGData: from :func:`SpcGroup`
1093    :param All: True return all equivalent positions including duplicates;
1094      False return only unique positions
1095    :param Uij: [U11,U22,U33,U12,U13,U23] or [] if no Uij
1096    :param Move: True move generated atom positions to be inside cell
1097      False do not move atoms       
1098    :return: [[XYZEquiv],Idup,[UijEquiv]]
1099
1100      *  [XYZEquiv] is list of equivalent positions (XYZ is first entry)
1101      *  Idup = [-][C]SS where SS is the symmetry operator number (1-24), C (if not 0,0,0)
1102      * is centering operator number (1-4) and - is for inversion
1103        Cell = unit cell translations needed to put new positions inside cell
1104        [UijEquiv] - equivalent Uij; absent if no Uij given
1105       
1106    '''
1107    XYZEquiv = []
1108    UijEquiv = []
1109    Idup = []
1110    Cell = []
1111    X = np.array(XYZ)
1112    celli = np.zeros(3)
1113    if Move:
1114        X,celli = MoveToUnitCell(X)
1115    for ic,cen in enumerate(SGData['SGCen']):
1116        C = np.array(cen)
1117        for invers in range(int(SGData['SGInv']+1)):
1118            for io,[M,T] in enumerate(SGData['SGOps']):
1119                idup = ((io+1)+100*ic)*(1-2*invers)
1120                XT = np.inner(M,X)+T
1121                if len(Uij):
1122                    U = Uij2U(Uij)
1123                    U = np.inner(M,np.inner(U,M).T)
1124                    newUij = U2Uij(U)
1125                if invers:
1126                    XT = -XT
1127                XT += C
1128                cell = np.zeros(3)
1129                cellj = np.zeros(3)
1130                if Move:
1131                    newX,cellj = MoveToUnitCell(XT)
1132                else:
1133                    newX = XT
1134                cell += (celli+cellj)
1135                if All:
1136                    if np.allclose(newX,X,atol=0.0002):
1137                        idup = False
1138                else:
1139                    if True in [np.allclose(newX,oldX,atol=0.0002) for oldX in XYZEquiv]:
1140                        idup = False
1141                if All or idup:
1142                    XYZEquiv.append(newX)
1143                    Idup.append(idup)
1144                    Cell.append(cell)
1145                    if len(Uij):
1146                        UijEquiv.append(newUij)                   
1147    if len(Uij):
1148        return zip(XYZEquiv,UijEquiv,Idup,Cell)
1149    else:
1150        return zip(XYZEquiv,Idup,Cell)
1151
1152def GenHKLf(HKL,SGData):
1153    '''
1154    Uses old GSAS Fortran routine genhkl.for
1155
1156    :param HKL:  [h,k,l] must be integral values for genhkl.for to work
1157    :param SGData: space group data obtained from SpcGroup
1158    :returns: iabsnt,mulp,Uniq,phi
1159
1160     *   iabsnt = True if reflection is forbidden by symmetry
1161     *   mulp = reflection multiplicity including Friedel pairs
1162     *   Uniq = numpy array of equivalent hkl in descending order of h,k,l
1163     *   phi = phase offset for each equivalent h,k,l
1164
1165    '''
1166    hklf = list(HKL)+[0,]       #could be numpy array!
1167    Ops = SGData['SGOps']
1168    OpM = np.array([op[0] for op in Ops])
1169    OpT = np.array([op[1] for op in Ops])
1170    Inv = SGData['SGInv']
1171    Cen = np.array([cen for cen in SGData['SGCen']])
1172   
1173    Nuniq,Uniq,iabsnt,mulp = pyspg.genhklpy(hklf,len(Ops),OpM,OpT,SGData['SGInv'],len(Cen),Cen)
1174    h,k,l,f = Uniq
1175    Uniq=np.array(zip(h[:Nuniq],k[:Nuniq],l[:Nuniq]))
1176    phi = f[:Nuniq]
1177   
1178    return iabsnt,mulp,Uniq,phi
1179   
1180def checkSSLaue(HKL,SGData,SSGData):
1181    #Laue check here - Toss HKL if outside unique Laue part
1182    h,k,l,m = HKL
1183    if SGData['SGLaue'] == '2/m':
1184        if SGData['SGUniq'] == 'a':
1185            if 'a' in SSGData['modSymb'] and h == 0 and m < 0:
1186                return False
1187            elif 'b' in SSGData['modSymb'] and k == 0 and l ==0 and m < 0:
1188                return False
1189            else:
1190                return True
1191        elif SGData['SGUniq'] == 'b':
1192            if 'b' in SSGData['modSymb'] and k == 0 and m < 0:
1193                return False
1194            elif 'a' in SSGData['modSymb'] and h == 0 and l ==0 and m < 0:
1195                return False
1196            else:
1197                return True
1198        elif SGData['SGUniq'] == 'c':
1199            if 'g' in SSGData['modSymb'] and l == 0 and m < 0:
1200                return False
1201            elif 'a' in SSGData['modSymb'] and h == 0 and k ==0 and m < 0:
1202                return False
1203            else:
1204                return True
1205    elif SGData['SGLaue'] == 'mmm':
1206        if 'a' in SSGData['modSymb']:
1207            if h == 0 and m < 0:
1208                return False
1209            else:
1210                return True
1211        elif 'b' in SSGData['modSymb']:
1212            if k == 0 and m < 0:
1213                return False
1214            else:
1215                return True
1216        elif 'g' in SSGData['modSymb']:
1217            if l == 0 and m < 0:
1218                return False
1219            else:
1220                return True
1221    else:   #tetragonal, trigonal, hexagonal (& triclinic?)
1222        if l == 0 and m < 0:
1223            return False
1224        else:
1225            return True
1226       
1227   
1228def checkSSextc(HKL,SSGData):
1229    Ops = SSGData['SSGOps']
1230    OpM = np.array([op[0] for op in Ops])
1231    OpT = np.array([op[1] for op in Ops])
1232    HKLS = np.array([HKL,-HKL])     #Freidel's Law
1233    DHKL = np.reshape(np.inner(HKLS,OpM)-HKL,(-1,4))
1234    PHKL = np.reshape(np.inner(HKLS,OpT),(-1,))
1235    for dhkl,phkl in zip(DHKL,PHKL)[1:]:    #skip identity
1236        if dhkl.any():
1237            continue
1238        else:
1239            if phkl%1.:
1240                return False
1241    return True
1242                                 
1243def GetOprPtrName(key):
1244    'Needs a doc string'
1245    OprPtrName = {
1246        '-6643':[   2,' 1bar ', 1],'6479' :[  10,'  2z  ', 2],'-6479':[   9,'  mz  ', 3],
1247        '6481' :[   7,'  my  ', 4],'-6481':[   6,'  2y  ', 5],'6641' :[   4,'  mx  ', 6],
1248        '-6641':[   3,'  2x  ', 7],'6591' :[  28,' m+-0 ', 8],'-6591':[  27,' 2+-0 ', 9],
1249        '6531' :[  25,' m110 ',10],'-6531':[  24,' 2110 ',11],'6537' :[  61,'  4z  ',12],
1250        '-6537':[  62,' -4z  ',13],'975'  :[  68,' 3+++1',14],'6456' :[ 114,'  3z1 ',15],
1251        '-489' :[  73,' 3+-- ',16],'483'  :[  78,' 3-+- ',17],'-969' :[  83,' 3--+ ',18],
1252        '819'  :[  22,' m+0- ',19],'-819' :[  21,' 2+0- ',20],'2431' :[  16,' m0+- ',21],
1253        '-2431':[  15,' 20+- ',22],'-657' :[  19,' m101 ',23],'657'  :[  18,' 2101 ',24],
1254        '1943' :[  48,' -4x  ',25],'-1943':[  47,'  4x  ',26],'-2429':[  13,' m011 ',27],
1255        '2429' :[  12,' 2011 ',28],'639'  :[  55,' -4y  ',29],'-639' :[  54,'  4y  ',30],
1256        '-6484':[ 146,' 2010 ', 4],'6484' :[ 139,' m010 ', 5],'-6668':[ 145,' 2100 ', 6],
1257        '6668' :[ 138,' m100 ', 7],'-6454':[ 148,' 2120 ',18],'6454' :[ 141,' m120 ',19],
1258        '-6638':[ 149,' 2210 ',20],'6638' :[ 142,' m210 ',21],              #search ends here
1259        '2223' :[  68,' 3+++2',39],
1260        '6538' :[ 106,'  6z1 ',40],'-2169':[  83,' 3--+2',41],'2151' :[  73,' 3+--2',42],
1261        '2205' :[  79,'-3-+-2',43],'-2205':[  78,' 3-+-2',44],'489'  :[  74,'-3+--1',45],
1262        '801'  :[  53,'  4y1 ',46],'1945' :[  47,'  4x3 ',47],'-6585':[  62,' -4z3 ',48],
1263        '6585' :[  61,'  4z3 ',49],'6584' :[ 114,'  3z2 ',50],'6666' :[ 106,'  6z5 ',51],
1264        '6643' :[   1,' Iden ',52],'-801' :[  55,' -4y1 ',53],'-1945':[  48,' -4x3 ',54],
1265        '-6666':[ 105,' -6z5 ',55],'-6538':[ 105,' -6z1 ',56],'-2223':[  69,'-3+++2',57],
1266        '-975' :[  69,'-3+++1',58],'-6456':[ 113,' -3z1 ',59],'-483' :[  79,'-3-+-1',60],
1267        '969'  :[  84,'-3--+1',61],'-6584':[ 113,' -3z2 ',62],'2169' :[  84,'-3--+2',63],
1268        '-2151':[  74,'-3+--2',64],'0':[0,' ????',0]
1269        }
1270    return OprPtrName[key]
1271
1272def GetKNsym(key):
1273    'Needs a doc string'
1274    KNsym = {
1275        '0'         :'    1   ','1'         :'   -1   ','64'        :'    2(x)','32'        :'    m(x)',
1276        '97'        :'  2/m(x)','16'        :'    2(y)','8'         :'    m(y)','25'        :'  2/m(y)',
1277        '2'         :'    2(z)','4'         :'    m(z)','7'         :'  2/m(z)','134217728' :'   2(yz)',
1278        '67108864'  :'   m(yz)','201326593' :' 2/m(yz)','2097152'   :'  2(0+-)','1048576'   :'  m(0+-)',
1279        '3145729'   :'2/m(0+-)','8388608'   :'   2(xz)','4194304'   :'   m(xz)','12582913'  :' 2/m(xz)',
1280        '524288'    :'  2(+0-)','262144'    :'  m(+0-)','796433'    :'2/m(+0-)','1024'      :'   2(xy)',
1281        '512'       :'   m(xy)','1537'      :' 2/m(xy)','256'       :'  2(+-0)','128'       :'  m(+-0)',
1282        '385'       :'2/m(+-0)','76'        :'  mm2(x)','52'        :'  mm2(y)','42'        :'  mm2(z)',
1283        '135266336' :' mm2(yz)','69206048'  :'mm2(0+-)','8650760'   :' mm2(xz)','4718600'   :'mm2(+0-)',
1284        '1156'      :' mm2(xy)','772'       :'mm2(+-0)','82'        :'  222   ','136314944' :'  222(x)',
1285        '8912912'   :'  222(y)','1282'      :'  222(z)','127'       :'  mmm   ','204472417' :'  mmm(x)',
1286        '13369369'  :'  mmm(y)','1927'      :'  mmm(z)','33554496'  :'  4(100)','16777280'  :' -4(100)',
1287        '50331745'  :'4/m(100)','169869394' :'422(100)','84934738'  :'-42m 100','101711948' :'4mm(100)',
1288        '254804095' :'4/mmm100','536870928 ':'  4(010)','268435472' :' -4(010)','805306393' :'4/m (10)',
1289        '545783890' :'422(010)','272891986' :'-42m 010','541327412' :'4mm(010)','818675839' :'4/mmm010',
1290        '2050'      :'  4(001)','4098'      :' -4(001)','6151'      :'4/m(001)','3410'      :'422(001)',
1291        '4818'      :'-42m 001','2730'      :'4mm(001)','8191'      :'4/mmm001','8192'      :'  3(111)',
1292        '8193'      :' -3(111)','2629888'   :' 32(111)','1319040'   :' 3m(111)','3940737'   :'-3m(111)',
1293        '32768'     :'  3(+--)','32769'     :' -3(+--)','10519552'  :' 32(+--)','5276160'   :' 3m(+--)',
1294        '15762945'  :'-3m(+--)','65536'     :'  3(-+-)','65537'     :' -3(-+-)','134808576' :' 32(-+-)',
1295        '67437056'  :' 3m(-+-)','202180097' :'-3m(-+-)','131072'    :'  3(--+)','131073'    :' -3(--+)',
1296        '142737664' :' 32(--+)','71434368'  :' 3m(--+)','214040961' :'-3m(--+)','237650'    :'   23   ',
1297        '237695'    :'   m3   ','715894098' :'   432  ','358068946' :'  -43m  ','1073725439':'   m3m  ',
1298        '68157504'  :' mm2d100','4456464'   :' mm2d010','642'       :' mm2d001','153092172' :'-4m2 100',
1299        '277348404' :'-4m2 010','5418'      :'-4m2 001','1075726335':'  6/mmm ','1074414420':'-6m2 100',
1300        '1075070124':'-6m2 120','1075069650':'   6mm  ','1074414890':'   622  ','1073758215':'   6/m  ',
1301        '1073758212':'   -6   ','1073758210':'    6   ','1073759865':'-3m(100)','1075724673':'-3m(120)',
1302        '1073758800':' 3m(100)','1075069056':' 3m(120)','1073759272':' 32(100)','1074413824':' 32(120)',
1303        '1073758209':'   -3   ','1073758208':'    3   ','1074135143':'mmm(100)','1075314719':'mmm(010)',
1304        '1073743751':'mmm(110)','1074004034':' mm2z100','1074790418':' mm2z010','1073742466':' mm2z110',
1305        '1074004004':'mm2(100)','1074790412':'mm2(010)','1073742980':'mm2(110)','1073872964':'mm2(120)',
1306        '1074266132':'mm2(210)','1073742596':'mm2(+-0)','1073872930':'222(100)','1074266122':'222(010)',
1307        '1073743106':'222(110)','1073741831':'2/m(001)','1073741921':'2/m(100)','1073741849':'2/m(010)',
1308        '1073743361':'2/m(110)','1074135041':'2/m(120)','1075314689':'2/m(210)','1073742209':'2/m(+-0)',
1309        '1073741828':' m(001) ','1073741888':' m(100) ','1073741840':' m(010) ','1073742336':' m(110) ',
1310        '1074003968':' m(120) ','1074790400':' m(210) ','1073741952':' m(+-0) ','1073741826':' 2(001) ',
1311        '1073741856':' 2(100) ','1073741832':' 2(010) ','1073742848':' 2(110) ','1073872896':' 2(120) ',
1312        '1074266112':' 2(210) ','1073742080':' 2(+-0) ','1073741825':'   -1   '
1313        }
1314    return KNsym[key]       
1315
1316def GetNXUPQsym(siteSym):
1317    '''       
1318    The codes XUPQ are for lookup of symmetry constraints for position(X), thermal parm(U) & magnetic moments
1319    (P&Q-not used in GSAS-II)
1320    '''
1321    NXUPQsym = {
1322        '    1   ':(28,29,28,28),'   -1   ':( 1,29,28, 0),'    2(x)':(12,18,12,25),'    m(x)':(25,18,12,25),
1323        '  2/m(x)':( 1,18, 0,-1),'    2(y)':(13,17,13,24),'    m(y)':(24,17,13,24),'  2/m(y)':( 1,17, 0,-1),
1324        '    2(z)':(14,16,14,23),'    m(z)':(23,16,14,23),'  2/m(z)':( 1,16, 0,-1),'   2(yz)':(10,23,10,22),
1325        '   m(yz)':(22,23,10,22),' 2/m(yz)':( 1,23, 0,-1),'  2(0+-)':(11,24,11,21),'  m(0+-)':(21,24,11,21),
1326        '2/m(0+-)':( 1,24, 0,-1),'   2(xz)':( 8,21, 8,20),'   m(xz)':(20,21, 8,20),' 2/m(xz)':( 1,21, 0,-1),
1327        '  2(+0-)':( 9,22, 9,19),'  m(+0-)':(19,22, 9,19),'2/m(+0-)':( 1,22, 0,-1),'   2(xy)':( 6,19, 6,18),
1328        '   m(xy)':(18,19, 6,18),' 2/m(xy)':( 1,19, 0,-1),'  2(+-0)':( 7,20, 7,17),'  m(+-0)':(17,20, 7,17),
1329        '2/m(+-0)':( 1,20, 0,-1),'  mm2(x)':(12,10, 0,-1),'  mm2(y)':(13,10, 0,-1),'  mm2(z)':(14,10, 0,-1),
1330        ' mm2(yz)':(10,13, 0,-1),'mm2(0+-)':(11,13, 0,-1),' mm2(xz)':( 8,12, 0,-1),'mm2(+0-)':( 9,12, 0,-1),
1331        ' mm2(xy)':( 6,11, 0,-1),'mm2(+-0)':( 7,11, 0,-1),'  222   ':( 1,10, 0,-1),'  222(x)':( 1,13, 0,-1),
1332        '  222(y)':( 1,12, 0,-1),'  222(z)':( 1,11, 0,-1),'  mmm   ':( 1,10, 0,-1),'  mmm(x)':( 1,13, 0,-1),
1333        '  mmm(y)':( 1,12, 0,-1),'  mmm(z)':( 1,11, 0,-1),'  4(100)':(12, 4,12, 0),' -4(100)':( 1, 4,12, 0),
1334        '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),
1335        '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),
1336        '422(010)':( 1, 3, 0,-1),'-42m 010':( 1, 3, 0,-1),'4mm(010)':(13, 3, 0,-1),'4/mmm010':(1, 3, 0,-1,),
1337        '  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),
1338        '-42m 001':( 1, 2, 0,-1),'4mm(001)':(14, 2, 0,-1),'4/mmm001':( 1, 2, 0,-1),'  3(111)':( 2, 5, 2, 0),
1339        ' -3(111)':( 1, 5, 2, 0),' 32(111)':( 1, 5, 0, 2),' 3m(111)':( 2, 5, 0, 2),'-3m(111)':( 1, 5, 0,-1),
1340        '  3(+--)':( 5, 8, 5, 0),' -3(+--)':( 1, 8, 5, 0),' 32(+--)':( 1, 8, 0, 5),' 3m(+--)':( 5, 8, 0, 5),
1341        '-3m(+--)':( 1, 8, 0,-1),'  3(-+-)':( 4, 7, 4, 0),' -3(-+-)':( 1, 7, 4, 0),' 32(-+-)':( 1, 7, 0, 4),
1342        ' 3m(-+-)':( 4, 7, 0, 4),'-3m(-+-)':( 1, 7, 0,-1),'  3(--+)':( 3, 6, 3, 0),' -3(--+)':( 1, 6, 3, 0),
1343        ' 32(--+)':( 1, 6, 0, 3),' 3m(--+)':( 3, 6, 0, 3),'-3m(--+)':( 1, 6, 0,-1),'   23   ':( 1, 1, 0, 0),
1344        '   m3   ':( 1, 1, 0, 0),'   432  ':( 1, 1, 0, 0),'  -43m  ':( 1, 1, 0, 0),'   m3m  ':( 1, 1, 0, 0),
1345        ' mm2d100':(12,13, 0,-1),' mm2d010':(13,12, 0,-1),' mm2d001':(14,11, 0,-1),'-4m2 100':( 1, 4, 0,-1),
1346        '-4m2 010':( 1, 3, 0,-1),'-4m2 001':( 1, 2, 0,-1),'  6/mmm ':( 1, 9, 0,-1),'-6m2 100':( 1, 9, 0,-1),
1347        '-6m2 120':( 1, 9, 0,-1),'   6mm  ':(14, 9, 0,-1),'   622  ':( 1, 9, 0,-1),'   6/m  ':( 1, 9,14,-1),
1348        '   -6   ':( 1, 9,14, 0),'    6   ':(14, 9,14, 0),'-3m(100)':( 1, 9, 0,-1),'-3m(120)':( 1, 9, 0,-1),
1349        ' 3m(100)':(14, 9, 0,14),' 3m(120)':(14, 9, 0,14),' 32(100)':( 1, 9, 0,14),' 32(120)':( 1, 9, 0,14),
1350        '   -3   ':( 1, 9,14, 0),'    3   ':(14, 9,14, 0),'mmm(100)':( 1,14, 0,-1),'mmm(010)':( 1,15, 0,-1),
1351        'mmm(110)':( 1,11, 0,-1),' mm2z100':(14,14, 0,-1),' mm2z010':(14,15, 0,-1),' mm2z110':(14,11, 0,-1),
1352        'mm2(100)':(12,14, 0,-1),'mm2(010)':(13,15, 0,-1),'mm2(110)':( 6,11, 0,-1),'mm2(120)':(15,14, 0,-1),
1353        'mm2(210)':(16,15, 0,-1),'mm2(+-0)':( 7,11, 0,-1),'222(100)':( 1,14, 0,-1),'222(010)':( 1,15, 0,-1),
1354        '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),
1355        '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),
1356        ' m(001) ':(23,16,14,23),' m(100) ':(26,25,12,26),' m(010) ':(27,28,13,27),' m(110) ':(18,19, 6,18),
1357        ' m(120) ':(24,27,15,24),' m(210) ':(25,26,16,25),' m(+-0) ':(17,20, 7,17),' 2(001) ':(14,16,14,23),
1358        ' 2(100) ':(12,25,12,26),' 2(010) ':(13,28,13,27),' 2(110) ':( 6,19, 6,18),' 2(120) ':(15,27,15,24),
1359        ' 2(210) ':(16,26,16,25),' 2(+-0) ':( 7,20, 7,17),'   -1   ':( 1,29,28, 0)
1360        }
1361    return NXUPQsym[siteSym]
1362
1363def GetCSxinel(siteSym): 
1364    'Needs a doc string'
1365    CSxinel = [[],                         # 0th empty - indices are Fortran style
1366        [[0,0,0],[ 0.0, 0.0, 0.0]],      #1  0  0  0
1367        [[1,1,1],[ 1.0, 1.0, 1.0]],      #2  X  X  X
1368        [[1,1,1],[ 1.0, 1.0,-1.0]],      #3  X  X -X
1369        [[1,1,1],[ 1.0,-1.0, 1.0]],      #4  X -X  X
1370        [[1,1,1],[ 1.0,-1.0,-1.0]],      #5 -X  X  X
1371        [[1,1,0],[ 1.0, 1.0, 0.0]],      #6  X  X  0
1372        [[1,1,0],[ 1.0,-1.0, 0.0]],      #7  X -X  0
1373        [[1,0,1],[ 1.0, 0.0, 1.0]],      #8  X  0  X
1374        [[1,0,1],[ 1.0, 0.0,-1.0]],      #9  X  0 -X
1375        [[0,1,1],[ 0.0, 1.0, 1.0]],      #10  0  Y  Y
1376        [[0,1,1],[ 0.0, 1.0,-1.0]],      #11 0  Y -Y
1377        [[1,0,0],[ 1.0, 0.0, 0.0]],      #12  X  0  0
1378        [[0,1,0],[ 0.0, 1.0, 0.0]],      #13  0  Y  0
1379        [[0,0,1],[ 0.0, 0.0, 1.0]],      #14  0  0  Z
1380        [[1,1,0],[ 1.0, 2.0, 0.0]],      #15  X 2X  0
1381        [[1,1,0],[ 2.0, 1.0, 0.0]],      #16 2X  X  0
1382        [[1,1,2],[ 1.0, 1.0, 1.0]],      #17  X  X  Z
1383        [[1,1,2],[ 1.0,-1.0, 1.0]],      #18  X -X  Z
1384        [[1,2,1],[ 1.0, 1.0, 1.0]],      #19  X  Y  X
1385        [[1,2,1],[ 1.0, 1.0,-1.0]],      #20  X  Y -X
1386        [[1,2,2],[ 1.0, 1.0, 1.0]],      #21  X  Y  Y
1387        [[1,2,2],[ 1.0, 1.0,-1.0]],      #22  X  Y -Y
1388        [[1,2,0],[ 1.0, 1.0, 0.0]],      #23  X  Y  0
1389        [[1,0,2],[ 1.0, 0.0, 1.0]],      #24  X  0  Z
1390        [[0,1,2],[ 0.0, 1.0, 1.0]],      #25  0  Y  Z
1391        [[1,1,2],[ 1.0, 2.0, 1.0]],      #26  X 2X  Z
1392        [[1,1,2],[ 2.0, 1.0, 1.0]],      #27 2X  X  Z
1393        [[1,2,3],[ 1.0, 1.0, 1.0]],      #28  X  Y  Z
1394        ]
1395    indx = GetNXUPQsym(siteSym)
1396    return CSxinel[indx[0]]
1397   
1398def GetCSuinel(siteSym):
1399    "returns Uij terms, multipliers, GUI flags & Uiso2Uij multipliers"
1400    CSuinel = [[],                                             # 0th empty - indices are Fortran style
1401        [[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]],    #1  A  A  A  0  0  0
1402        [[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]],    #2  A  A  C  0  0  0
1403        [[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]],    #3  A  B  A  0  0  0
1404        [[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]],    #4  A  B  B  0  0  0
1405        [[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]],    #5  A  A  A  D  D  D
1406        [[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]],    #6  A  A  A  D -D -D
1407        [[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]],    #7  A  A  A  D -D  D
1408        [[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]],    #8  A  A  A  D  D -D
1409        [[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]],    #9  A  A  C A/2 0  0
1410        [[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]],    #10  A  B  C  0  0  0
1411        [[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]],    #11  A  A  C  D  0  0
1412        [[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]],    #12  A  B  A  0  E  0
1413        [[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]],    #13  A  B  B  0  0  F
1414        [[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]],    #14  A  B  C B/2 0  0
1415        [[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]],    #15  A  B  C A/2 0  0
1416        [[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]],    #16  A  B  C  D  0  0
1417        [[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]],    #17  A  B  C  0  E  0
1418        [[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]],    #18  A  B  C  0  0  F
1419        [[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]],    #19  A  A  C  D  E -E
1420        [[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]],    #20  A  A  C  D  E  E
1421        [[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]],    #21  A  B  A  D  E -D
1422        [[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]],    #22  A  B  A  D  E  D
1423        [[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]],    #23  A  B  B  D -D  F
1424        [[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]],    #24  A  B  B  D  D  F
1425        [[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]],    #25  A  B  C B/2 F/2 F
1426        [[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]],    #26  A  B  C A/2  0  F
1427        [[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]],    #27  A  B  C B/2  E  0
1428        [[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]],    #28  A  B  C A/2  E E/2
1429        [[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]],    #29  A  B  C  D  E   F
1430        ]
1431    indx = GetNXUPQsym(siteSym)
1432    return CSuinel[indx[1]]
1433   
1434#def GetSSfxuinel(waveType,nH,XYZ,SGData,SSGData,debug=False):
1435#   
1436#    def fracCrenel(tau,Toff,Twid):
1437#        Tau = (tau-Toff)%1.
1438#        A = np.where(Tau<Twid,1.,0.)
1439#        return A
1440#       
1441#    def fracFourier(tau,nH,fsin,fcos):
1442#        SA = np.sin(2.*nH*np.pi*tau)
1443#        CB = np.cos(2.*nH*np.pi*tau)
1444#        A = SA[np.newaxis,np.newaxis,:]*fsin[:,:,np.newaxis]
1445#        B = CB[np.newaxis,np.newaxis,:]*fcos[:,:,np.newaxis]
1446#        return A+B
1447#       
1448#    def posFourier(tau,nH,psin,pcos):
1449#        SA = np.sin(2*nH*np.pi*tau)
1450#        CB = np.cos(2*nH*np.pi*tau)
1451#        A = SA[np.newaxis,np.newaxis,:]*psin[:,:,np.newaxis]
1452#        B = CB[np.newaxis,np.newaxis,:]*pcos[:,:,np.newaxis]
1453#        return A+B   
1454#
1455#    def posSawtooth(tau,Toff,slopes):
1456#        Tau = (tau-Toff[:,np.newaxis])%1.
1457#        A = slopes[:,:,np.newaxis]*Tau
1458#        return A
1459#   
1460#    def posZigZag(tau,Toff,slopes):
1461#        Tau = (tau-Toff[:,np.newaxis])%1.
1462#        A = np.where(Tau <= 0.5,slopes[:,:,np.newaxis]*Tau,slopes[:,:,np.newaxis]*(1.-Tau))
1463#        return A
1464#       
1465#    print 'super space group: ',SSGData['SSpGrp']
1466#    CSI = {'Sfrac':[[[1,0],[2,0]],[[1.,0.],[1.,0.]]],
1467#        'Spos':[[[1,0,0],[2,0,0],[3,0,0], [4,0,0],[5,0,0],[6,0,0]],
1468#            [[1.,0.,0.],[1.,0.,0.],[1.,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.]]],    #sin & cos
1469#        'Sadp':[[[1,0,0],[2,0,0],[3,0,0],[4,0,0],[5,0,0],[6,0,0],
1470#            [7,0,0],[8,0,0],[9,0,0],[10,0,0],[11,0,0],[12,0,0]],
1471#            [[1.,0.,0.],[1.,0.,0.],[1.,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.],
1472#            [1.,0.,0.],[1.,0.,0.],[1.,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.]]],
1473#        'Smag':[[[1,0,0],[2,0,0],[3,0,0], [4,0,0],[5,0,0],[6,0,0]],
1474#            [[1.,0.,0.],[1.,0.,0.],[1.,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.]]],}
1475#    xyz = np.array(XYZ)%1.
1476#    xyzt = np.array(XYZ+[0,])%1.
1477#    SGOps = copy.deepcopy(SGData['SGOps'])
1478#    siteSym = SytSym(XYZ,SGData)[0].strip()
1479#    print 'siteSym: ',siteSym
1480#    if siteSym == '1':   #"1" site symmetry
1481#        if debug:
1482#            return CSI,None,None,None,None
1483#        else:
1484#            return CSI
1485#    elif siteSym == '-1':   #"-1" site symmetry
1486#        CSI['Sfrac'][0] = [[1,0],[0,0]]
1487#        CSI['Spos'][0] = [[1,0,0],[2,0,0],[3,0,0], [0,0,0],[0,0,0],[0,0,0]]
1488#        CSI['Sadp'][0] = [[0,0,0],[0,0,0],[0,0,0],[0,0,0],[0,0,0],[0,0,0],
1489#        [1,0,0],[2,0,0],[3,0,0],[4,0,0],[5,0,0],[6,0,0]]
1490#        if debug:
1491#            return CSI,None,None,None,None
1492#        else:
1493#            return CSI
1494#    SSGOps = copy.deepcopy(SSGData['SSGOps'])
1495#    #expand ops to include inversions if any
1496#    if SGData['SGInv']:
1497#        for op,sop in zip(SGData['SGOps'],SSGData['SSGOps']):
1498#            SGOps.append([-op[0],-op[1]%1.])
1499#            SSGOps.append([-sop[0],-sop[1]%1.])
1500#    #build set of sym ops around special poasition       
1501#    SSop = []
1502#    Sop = []
1503#    Sdtau = []
1504#    for iop,Op in enumerate(SGOps):         
1505#        nxyz = (np.inner(Op[0],xyz)+Op[1])%1.
1506#        if np.allclose(xyz,nxyz,1.e-4) and iop and MT2text(Op).replace(' ','') != '-X,-Y,-Z':
1507#            SSop.append(SSGOps[iop])
1508#            Sop.append(SGOps[iop])
1509#            ssopinv = nl.inv(SSGOps[iop][0])
1510#            mst = ssopinv[3][:3]
1511#            epsinv = ssopinv[3][3]
1512#            Sdtau.append(np.sum(mst*(XYZ-SGOps[iop][1])-epsinv*SSGOps[iop][1][3]))
1513#    SdIndx = np.argsort(np.array(Sdtau))     # just to do in sensible order
1514#    OpText =  [MT2text(s).replace(' ','') for s in Sop]         #debug?
1515#    SSOpText = [SSMT2text(ss).replace(' ','') for ss in SSop]   #debug?
1516#    print 'special pos super operators: ',SSOpText
1517#    #setup displacement arrays
1518#    tau = np.linspace(0,1,49,True)
1519#    delt2 = np.eye(2)*0.001
1520#    delt4 = np.eye(4)*0.001
1521#    delt6 = np.eye(6)*0.001
1522#    delt12 = np.eye(12)*0.0001
1523#    #make modulation arrays - one parameter at a time
1524#    #site fractions
1525#    CSI['Sfrac'] = [np.zeros((2),dtype='i'),np.ones(2)]
1526#    if 'Crenel' in waveType:
1527#        dF = fracCrenel(tau,delt2[:1],delt2[1:]).squeeze()
1528#    else:
1529#        dF = fracFourier(tau,nH,delt2[:1],delt2[1:]).squeeze()
1530#    dFT = np.zeros_like(dF)
1531#    #positions       
1532#    if 'Fourier' in waveType:
1533#        dX = posFourier(tau,nH,delt6[:3],delt6[3:]) #+np.array(XYZ)[:,np.newaxis,np.newaxis]
1534#          #3x6x12 modulated position array (X,Spos,tau)& force positive
1535#        CSI['Spos'] = [np.zeros((6,3),dtype='i'),np.zeros((6,3))]
1536#    elif waveType == 'Sawtooth':
1537#        CSI['Spos'] = [np.array([[1,],[2,],[3,],[4,]]),np.array([[1.0,],[1.0,],[1.0,],[1.0,]])]
1538#    elif waveType == 'ZigZag':
1539#        CSI['Spos'] = [np.array([[1,],[2,],[3,],[4,]]),np.array([[1.0,],[1.0,],[1.0,],[1.0,]])]
1540#    #anisotropic thermal motion
1541#    dU = posFourier(tau,nH,delt12[:6],delt12[6:])                  #Uij modulations - 6x12x12 array
1542#    CSI['Sadp'] = [np.zeros((12,3),dtype='i'),np.zeros((12,3))]
1543#       
1544#    FSC = np.ones(2,dtype='i')
1545#    VFSC = np.ones(2)
1546#    XSC = np.ones(6,dtype='i')
1547#    USC = np.ones(12,dtype='i')
1548#    dFTP = []
1549#    dXTP = []
1550#    dUTP = []
1551#    for i in SdIndx:
1552#        sop = Sop[i]
1553#        ssop = SSop[i]
1554#        fsc = np.ones(2,dtype='i')
1555#        xsc = np.ones(6,dtype='i')
1556#        ssopinv = nl.inv(ssop[0])
1557#        mst = ssopinv[3][:3]
1558#        epsinv = ssopinv[3][3]
1559#        sdet = nl.det(sop[0])
1560#        ssdet = nl.det(ssop[0])
1561#        dtau = mst*(XYZ-sop[1])-epsinv*ssop[1][3]
1562#        dT = 1.0
1563#        if np.any(dtau%.5):
1564#            dT = np.tan(np.pi*np.sum(dtau%.5))
1565#        tauT = np.inner(mst,XYZ-sop[1])+epsinv*(tau-ssop[1][3])
1566#        if waveType == 'Fourier':
1567#            dXT = posFourier(np.sort(tauT),nH,delt6[:3],delt6[3:])   #+np.array(XYZ)[:,np.newaxis,np.newaxis]
1568#        elif waveType == 'Sawtooth':
1569#            dXT = posSawtooth(tauT,delt4[0],delt4[1:])+np.array(XYZ)[:,np.newaxis,np.newaxis]
1570#        elif waveType == 'ZigZag':
1571#            dXT = posZigZag(tauT,delt4[0],delt4[1:])+np.array(XYZ)[:,np.newaxis,np.newaxis]           
1572#        dXT = np.inner(sop[0],dXT.T)
1573#        dXT = np.swapaxes(dXT,1,2)
1574#        dXT[:,:3,:] *= ssdet
1575#        dXTP.append(dXT)
1576#        if waveType == 'Fourier':
1577#            if np.any(dtau%.5) and ('1/2' in SSGData['modSymb'] or '1' in SSGData['modSymb']):
1578#                CSI['Spos'] = [[[1,0,0],[2,0,0],[3,0,0], [1,0,0],[2,0,0],[3,0,0]],
1579#                    [[1.,0.,0.],[1.,0.,0.],[1.,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.]]]                   
1580#                if '(x)' in siteSym:
1581#                    CSI['Spos'][1][3:] = [1./dT,0.,0.],[-dT,0.,0.],[-dT,0.,0.]
1582#                    if 'm' in siteSym and len(SdIndx) == 1:
1583#                        CSI['Spos'][1][3:] = [-dT,0.,0.],[1./dT,0.,0.],[1./dT,0.,0.]
1584#                elif '(y)' in siteSym:
1585#                    CSI['Spos'][1][3:] = [-dT,0.,0.],[1./dT,0.,0.],[-dT,0.,0.]
1586#                    if 'm' in siteSym and len(SdIndx) == 1:
1587#                        CSI['Spos'][1][3:] = [1./dT,0.,0.],[-dT,0.,0.],[1./dT,0.,0.]
1588#                elif '(z)' in siteSym:
1589#                    CSI['Spos'][1][3:] = [-dT,0.,0.],[-dT,0.,0.],[1./dT,0.,0.]
1590#                    if 'm' in siteSym and len(SdIndx) == 1:
1591#                        CSI['Spos'][1][3:] = [1./dT,0.,0.],[1./dT,0.,0.],[-dT,0.,0.]
1592#                for i in range(3):
1593#                    if not XSC[i]:
1594#                        CSI['Spos'][0][i] = [0,0,0]
1595#                        CSI['Spos'][1][i] = [0.,0.,0.]
1596#                        CSI['Spos'][0][i+3] = [0,0,0]
1597#                        CSI['Spos'][1][i+3] = [0.,0.,0.]
1598#            else:
1599#                for i in range(3):
1600#                    if np.allclose(dX[i,i,:],dXT[i,i,:]*sdet):
1601#                        xsc[i] = 1
1602#                    else:
1603#                        xsc[i] = 0
1604#                    if np.allclose(dX[i,i+3,:],dXT[i,i+3,:]):
1605#                        xsc[i+3] = 1
1606#                    else:
1607#                        xsc[i+3] = 0
1608#            XSC &= xsc
1609#           
1610#        fsc = np.ones(2,dtype='i')
1611#        if 'Crenel' in waveType:
1612#            dFT = fracCrenel(tauT,delt2[:1],delt2[1:]).squeeze()
1613#            fsc = [1,1]
1614#        else:
1615#            dFT = fracFourier(tauT,nH,delt2[:1],delt2[1:]).squeeze()
1616#            dFT = nl.det(sop[0])*dFT
1617#            dFT = dFT[:,np.argsort(tauT)]
1618#            dFT[0] *= ssdet
1619#            dFT[1] *= sdet
1620#            dFTP.append(dFT)
1621#       
1622#            if np.any(dtau%.5) and ('1/2' in SSGData['modSymb'] or '1' in SSGData['modSymb']):
1623#                fsc = [1,1]
1624#                CSI['Sfrac'] = [[[1,0],[1,0]],[[1.,0.],[1/dT,0.]]]
1625#                for i in range(2):
1626#                    if not FSC[i]:
1627#                        CSI['Sfrac'][0][i] = [0,0]
1628#                        CSI['Sfrac'][1][i] = [0.,0.]
1629#            else:
1630#                for i in range(2):
1631#                    if np.allclose(dF[i,:],dFT[i,:],atol=1.e-6):
1632#                        fsc[i] = 1
1633#                    else:
1634#                        fsc[i] = 0
1635#        FSC &= fsc
1636#           
1637#        usc = np.ones(12,dtype='i')
1638#        dUT = posFourier(tauT,nH,delt12[:6],delt12[6:])                  #Uij modulations - 6x12x49 array
1639#        dUijT = np.rollaxis(np.rollaxis(np.array(Uij2U(dUT)),3),3)    #convert dUT to 12x49x3x3
1640#        dUijT = np.rollaxis(np.inner(np.inner(sop[0],dUijT),sop[0].T),3)
1641#        dUT = np.array(U2Uij(dUijT))
1642#        dUT = dUT[:,:,np.argsort(tauT)]
1643#        dUTP.append(dUT)
1644#        if np.any(dtau%.5) and ('1/2' in SSGData['modSymb'] or '1' in SSGData['modSymb']):
1645#            CSI['Sadp'] = [[[1,0,0],[2,0,0],[3,0,0],[4,0,0],[5,0,0],[6,0,0],
1646#            [1,0,0],[2,0,0],[3,0,0],[4,0,0],[5,0,0],[6,0,0]],
1647#            [[1.,0.,0.],[1.,0.,0.],[1.,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.],
1648#            [1./dT,0.,0.],[1./dT,0.,0.],[1./dT,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.]]]
1649#            if '(x)' in siteSym:
1650#                CSI['Sadp'][1][9:] = [-dT,0.,0.],[-dT,0.,0.],[1./dT,0.,0.]
1651#            elif '(y)' in siteSym:
1652#                CSI['Sadp'][1][9:] = [-dT,0.,0.],[1./dT,0.,0.],[-dT,0.,0.]
1653#            elif '(z)' in siteSym:
1654#                CSI['Sadp'][1][9:] = [1./dT,0.,0.],[-dT,0.,0.],[-dT,0.,0.]
1655#            for i in range(6):
1656#                if not USC[i]:
1657#                    CSI['Sadp'][0][i] = [0,0,0]
1658#                    CSI['Sadp'][1][i] = [0.,0.,0.]
1659#                    CSI['Sadp'][0][i+6] = [0,0,0]
1660#                    CSI['Sadp'][1][i+6] = [0.,0.,0.]
1661#        else:
1662#                       
1663#            for i in range(6):
1664#                if np.allclose(dU[i,i,:],dUT[i,i,:]*sdet):
1665#                    usc[i] = 1
1666#                else:
1667#                    usc[i] = 0
1668#                if np.allclose(dU[i,i+6,:],dUT[i,i+6,:]):
1669#                    usc[i+6] = 1
1670#                else:
1671#                    usc[i+6] = 0
1672#            if '4/m' in siteSym and np.any(dUT[0,1,:]):
1673#                CSI['Sadp'][0][6:8] = [[12,0,0],[12,0,0]]
1674#                if ssop[1][3]:
1675#                    CSI['Sadp'][1][6:8] = [[1.,0.,0.],[-1.,0.,0.]]
1676#                    usc[9] = 1
1677#                else:
1678#                    CSI['Sadp'][1][6:8] = [[1.,0.,0.],[1.,0.,0.]]
1679#                    usc[9] = 0
1680#            elif '4' in siteSym and np.any(dUT[0,1,:]):
1681#                CSI['Sadp'][0][6:8] = [[12,0,0],[12,0,0]]
1682#                CSI['Sadp'][0][:2] = [[11,0,0],[11,0,0]]
1683#                if ssop[1][3]:
1684#                    CSI['Sadp'][1][:2] = [[1.,0.,0.],[-1.,0.,0.]]
1685#                    CSI['Sadp'][1][6:8] = [[1.,0.,0.],[-1.,0.,0.]]
1686#                    usc[3] = 1
1687#                    usc[9] = 1
1688#                else:
1689#                    CSI['Sadp'][1][:2] = [[1.,0.,0.],[1.,0.,0.]]
1690#                    CSI['Sadp'][1][6:8] = [[1.,0.,0.],[1.,0.,0.]]
1691#                    usc[3] = 0               
1692#                    usc[9] = 0
1693#            print SSMT2text(ssop).replace(' ',''),sdet,ssdet,epsinv,usc
1694#        USC &= usc
1695#    print dtau
1696#    print FSC
1697#    print XSC
1698#    print USC
1699#    if not np.any(dtau%.5):
1700#        n = -1
1701#        for i,U in enumerate(USC):
1702#            if U:
1703#                n += 1
1704#                CSI['Sadp'][0][i][0] = n+1
1705#                CSI['Sadp'][1][i][0] = 1.0
1706#        if waveType == 'Fourier':
1707#            n = -1
1708#            for i,X in enumerate(XSC):
1709#                if X:
1710#                    n += 1
1711#                    CSI['Spos'][0][i][0] = n+1
1712#                    CSI['Spos'][1][i][0] = 1.0
1713#        n = -1
1714#        for i,F in enumerate(FSC):
1715#            if F:
1716#                n += 1
1717#                CSI['Sfrac'][0][i] = n+1
1718#                CSI['Sfrac'][1][i] = 1.0
1719#            else:
1720#                CSI['Sfrac'][0][i] = 0
1721#                CSI['Sfrac'][1][i] = 0.           
1722#    if debug:
1723#        return CSI,[tau,tauT],[dF,dFTP],[dX,dXTP],[dU,dUTP]
1724#    else:
1725#        return CSI
1726#
1727def GetSSfxuinel(waveType,nH,XYZ,SGData,SSGData,debug=False):
1728   
1729    def orderParms(CSI):
1730        parms = [0,]
1731        for csi in CSI:
1732            for i in [0,1,2]:
1733                if csi[i] not in parms:
1734                    parms.append(csi[i])
1735        for csi in CSI:
1736            for i in [0,1,2]:
1737                csi[i] = parms.index(csi[i])
1738        return CSI
1739   
1740    def fracCrenel(tau,Toff,Twid):
1741        Tau = (tau-Toff[:,np.newaxis])%1.
1742        A = np.where(Tau<Twid[:,np.newaxis],1.,0.)
1743        return A
1744       
1745    def fracFourier(tau,nH,fsin,fcos):
1746        SA = np.sin(2.*nH*np.pi*tau)
1747        CB = np.cos(2.*nH*np.pi*tau)
1748        A = SA[np.newaxis,np.newaxis,:]*fsin[:,:,np.newaxis]
1749        B = CB[np.newaxis,np.newaxis,:]*fcos[:,:,np.newaxis]
1750        return A+B
1751       
1752    def posFourier(tau,nH,psin,pcos):
1753        SA = np.sin(2*nH*np.pi*tau)
1754        CB = np.cos(2*nH*np.pi*tau)
1755        A = SA[np.newaxis,np.newaxis,:]*psin[:,:,np.newaxis]
1756        B = CB[np.newaxis,np.newaxis,:]*pcos[:,:,np.newaxis]
1757        return A+B   
1758
1759    def posSawtooth(tau,Toff,slopes):
1760        Tau = (tau-Toff)%1.
1761        A = slopes[:,np.newaxis]*Tau
1762        return A
1763   
1764    def posZigZag(tau,Toff,slopes):
1765        Tau = (tau-Toff)%1.
1766        A = np.where(Tau <= 0.5,slopes[:,np.newaxis]*Tau,slopes[:,np.newaxis]*(1.-Tau))
1767        return A
1768       
1769    def DoFrac():
1770        delt2 = np.eye(2)*0.001
1771        FSC = np.ones(2,dtype='i')
1772        VFSC = np.ones(2)
1773        CSI = [np.zeros((2),dtype='i'),np.zeros(2)]
1774        if 'Crenel' in waveType:
1775            dF = fracCrenel(tau,delt2[:1],delt2[1:]).squeeze()
1776        else:
1777            dF = fracFourier(tau,nH,delt2[:1],delt2[1:]).squeeze()
1778        dFT = np.zeros_like(dF)
1779        dFTP = []
1780        for i in SdIndx:
1781            sop = Sop[i]
1782            ssop = SSop[i]
1783            ssopinv = nl.inv(ssop[0])
1784            mst = ssopinv[3][:3]
1785            epsinv = ssopinv[3][3]
1786            sdet = nl.det(sop[0])
1787            ssdet = nl.det(ssop[0])
1788            dtau = mst*(XYZ-sop[1])-epsinv*ssop[1][3]
1789            dT = 1.0
1790            if np.any(dtau%.5):
1791                dT = np.tan(np.pi*np.sum(dtau%.5))
1792            tauT = np.inner(mst,XYZ-sop[1])+epsinv*(tau-ssop[1][3])
1793            fsc = np.ones(2,dtype='i')
1794            if 'Crenel' in waveType:
1795                dFT = fracCrenel(tauT,delt2[:1],delt2[1:]).squeeze()
1796                fsc = [1,1]
1797            else:   #Fourier
1798                dFT = fracFourier(tauT,nH,delt2[:1],delt2[1:]).squeeze()
1799                dFT = nl.det(sop[0])*dFT
1800                dFT = dFT[:,np.argsort(tauT)]
1801                dFT[0] *= ssdet
1802                dFT[1] *= sdet
1803                dFTP.append(dFT)
1804           
1805                if np.any(dtau%.5) and ('1/2' in SSGData['modSymb'] or '1' in SSGData['modSymb']):
1806                    fsc = [1,1]
1807                    CSI = [[[1,0],[1,0]],[[1.,0.],[1/dT,0.]]]
1808                    FSC = np.zeros(2,dtype='i')
1809                    return CSI,dF,dFTP
1810                else:
1811                    for i in range(2):
1812                        if np.allclose(dF[i,:],dFT[i,:],atol=1.e-6):
1813                            fsc[i] = 1
1814                        else:
1815                            fsc[i] = 0
1816                    FSC &= fsc
1817                    print SSMT2text(ssop).replace(' ',''),sdet,ssdet,epsinv,fsc
1818        n = -1
1819        for i,F in enumerate(FSC):
1820            if F:
1821                n += 1
1822                CSI[0][i] = n+1
1823                CSI[1][i] = 1.0
1824       
1825        return CSI,dF,dFTP
1826       
1827    def DoXYZ():
1828        delt4 = np.ones(4)*0.001
1829        delt6 = np.eye(6)*0.001
1830        if 'Fourier' in waveType:
1831            dX = posFourier(tau,nH,delt6[:3],delt6[3:]) #+np.array(XYZ)[:,np.newaxis,np.newaxis]
1832              #3x6x12 modulated position array (X,Spos,tau)& force positive
1833            CSI = [np.zeros((6,3),dtype='i'),np.zeros((6,3))]
1834        elif waveType == 'Sawtooth':
1835            dX = posSawtooth(tau,delt4[0],delt4[1:])
1836            CSI = [np.array([[1,0,0],[2,0,0],[3,0,0],[4,0,0]]),
1837                np.array([[1.0,.0,.0],[1.0,.0,.0],[1.0,.0,.0],[1.0,.0,.0]])]
1838        elif waveType == 'ZigZag':
1839            dX = posZigZag(tau,delt4[0],delt4[1:])
1840            CSI = [np.array([[1,0,0],[2,0,0],[3,0,0],[4,0,0]]),
1841                np.array([[1.0,.0,.0],[1.0,.0,.0],[1.0,.0,.0],[1.0,.0,.0]])]
1842        XSC = np.ones(6,dtype='i')
1843        dXTP = []
1844        for i in SdIndx:
1845            sop = Sop[i]
1846            ssop = SSop[i]
1847            xsc = np.ones(6,dtype='i')
1848            ssopinv = nl.inv(ssop[0])
1849            mst = ssopinv[3][:3]
1850            epsinv = ssopinv[3][3]
1851            sdet = nl.det(sop[0])
1852            ssdet = nl.det(ssop[0])
1853            dtau = mst*(XYZ-sop[1])-epsinv*ssop[1][3]
1854            dT = 1.0
1855            if np.any(dtau%.5):
1856                dT = np.tan(np.pi*np.sum(dtau%.5))
1857            tauT = np.inner(mst,XYZ-sop[1])+epsinv*(tau-ssop[1][3])
1858            if waveType == 'Fourier':
1859                dXT = posFourier(np.sort(tauT),nH,delt6[:3],delt6[3:])   #+np.array(XYZ)[:,np.newaxis,np.newaxis]
1860            elif waveType == 'Sawtooth':
1861                dXT = posSawtooth(tauT,delt4[0],delt4[1:])+np.array(XYZ)[:,np.newaxis,np.newaxis]
1862            elif waveType == 'ZigZag':
1863                dXT = posZigZag(tauT,delt4[0],delt4[1:])+np.array(XYZ)[:,np.newaxis,np.newaxis] 
1864            dXT = np.inner(sop[0],dXT.T)    # X modulations array(3x6x49) -> array(3x49x6)
1865            dXT = np.swapaxes(dXT,1,2)      # back to array(3x6x49)
1866            dXT[:,:3,:] *= (ssdet*sdet)            # modify the sin component
1867            dXTP.append(dXT)
1868#            print np.sum(dX**2,axis=-1)*1e5
1869#            print np.sum(dXT**2,axis=-1)*1e5
1870            if waveType == 'Fourier':
1871                if np.any(dtau%.5) and ('1/2' in SSGData['modSymb'] or '1' in SSGData['modSymb']):
1872                    xsc[3:6] = 0
1873                    CSI = [[[1,0,0],[2,0,0],[3,0,0], [1,0,0],[2,0,0],[3,0,0]],
1874                        [[1.,0.,0.],[1.,0.,0.],[1.,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.]]]                   
1875                    if '(x)' in siteSym:
1876                        CSI[1][3:] = [1./dT,0.,0.],[-dT,0.,0.],[-dT,0.,0.]
1877                        if 'm' in siteSym and len(SdIndx) == 1:
1878                            CSI[1][3:] = [-dT,0.,0.],[1./dT,0.,0.],[1./dT,0.,0.]
1879                    elif '(y)' in siteSym:
1880                        CSI[1][3:] = [-dT,0.,0.],[1./dT,0.,0.],[-dT,0.,0.]
1881                        if 'm' in siteSym and len(SdIndx) == 1:
1882                            CSI[1][3:] = [1./dT,0.,0.],[-dT,0.,0.],[1./dT,0.,0.]
1883                    elif '(z)' in siteSym:
1884                        CSI[1][3:] = [-dT,0.,0.],[-dT,0.,0.],[1./dT,0.,0.]
1885                        if 'm' in siteSym and len(SdIndx) == 1:
1886                            CSI[1][3:] = [1./dT,0.,0.],[1./dT,0.,0.],[-dT,0.,0.]
1887                    elif '(xy)' in siteSym:
1888                        CSI[0] = [[1,0,0],[1,0,0],[2,0,0], [1,0,0],[1,0,0],[2,0,0]]
1889                        CSI[1][3:] = [[1./dT,0.,0.],[1./dT,0.,0.],[-dT,0.,0.]]
1890                elif '(xy)' in siteSym or '(+-0)' in siteSym:
1891                    mul = 1
1892                    if '(xy)' in siteSym:
1893                        mul = -1
1894                    if np.allclose(dX[0,0,:],dXT[1,0,:]*mul):
1895                        CSI[0][3:5] = [[11,0,0],[11,0,0]]
1896                        CSI[1][3:5] = [[1.,0,0],[mul,0,0]]
1897                        xsc[3:5] = 0
1898                    if np.allclose(dX[0,3,:],dXT[0,4,:]*mul):
1899                        CSI[0][:2] = [[12,0,0],[12,0,0]]
1900                        CSI[1][:2] = [[1.,0,0],[mul,0,0]]
1901                        xsc[:2] = 0
1902                else:
1903                    for i in range(3):
1904                        if not np.allclose(dX[i,i,:],dXT[i,i,:]):
1905                            xsc[i] = 0
1906                        if not np.allclose(dX[i,i+3,:],dXT[i,i+3,:]):
1907                            xsc[i+3] = 0
1908                XSC &= xsc
1909                print SSMT2text(ssop).replace(' ',''),sdet,ssdet,epsinv,xsc
1910        if waveType == 'Fourier':
1911            n = -1
1912            print XSC
1913            for i,X in enumerate(XSC):
1914                if X:
1915                    n += 1
1916                    CSI[0][i][0] = n+1
1917                    CSI[1][i][0] = 1.0
1918       
1919        return CSI,dX,dXTP
1920       
1921    def DoUij():
1922        tau = np.linspace(0,1,49,True)
1923        delt12 = np.eye(12)*0.0001
1924        dU = posFourier(tau,nH,delt12[:6],delt12[6:])                  #Uij modulations - 6x12x12 array
1925        CSI = [np.zeros((12,3),dtype='i'),np.zeros((12,3))]
1926        USC = np.ones(12,dtype='i')
1927        dUTP = []
1928        for i in SdIndx:
1929            sop = Sop[i]
1930            ssop = SSop[i]
1931            ssopinv = nl.inv(ssop[0])
1932            mst = ssopinv[3][:3]
1933            epsinv = ssopinv[3][3]
1934            sdet = nl.det(sop[0])
1935            ssdet = nl.det(ssop[0])
1936            dtau = mst*(XYZ-sop[1])-epsinv*ssop[1][3]
1937            dT = 1.0
1938            if np.any(dtau%.5):
1939                dT = np.tan(np.pi*np.sum(dtau%.5))
1940            tauT = np.inner(mst,XYZ-sop[1])+epsinv*(tau-ssop[1][3])
1941            usc = np.ones(12,dtype='i')
1942            dUT = posFourier(tauT,nH,delt12[:6],delt12[6:])                  #Uij modulations - 6x12x49 array
1943            dUijT = np.rollaxis(np.rollaxis(np.array(Uij2U(dUT)),3),3)    #convert dUT to 12x49x3x3
1944            dUijT = np.rollaxis(np.inner(np.inner(sop[0],dUijT),sop[0].T),3) #transform by sop - 3x3x12x49
1945            dUT = np.array(U2Uij(dUijT))    #convert to 6x12x49
1946            dUT = dUT[:,:,np.argsort(tauT)]
1947            dUTP.append(dUT)
1948            if np.any(dtau%.5) and ('1/2' in SSGData['modSymb'] or '1' in SSGData['modSymb']):
1949                CSI = [[[1,0,0],[2,0,0],[3,0,0],[4,0,0],[5,0,0],[6,0,0], 
1950                [1,0,0],[2,0,0],[3,0,0],[4,0,0],[5,0,0],[6,0,0]],
1951                [[1.,0.,0.],[1.,0.,0.],[1.,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.],
1952                [1./dT,0.,0.],[1./dT,0.,0.],[1./dT,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.]]]
1953                if 'mm2(x)' in siteSym:
1954                    CSI[1][9:] = [0.,0.,0.],[-dT,0.,0.],[0.,0.,0.]
1955                    USC = [1,1,1,0,1,0,1,1,1,0,1,0]
1956                elif '(xy)' in siteSym:
1957                    CSI[0] = [[1,0,0],[1,0,0],[2,0,0],[3,0,0],[4,0,0],[4,0,0],
1958                        [1,0,0],[1,0,0],[2,0,0],[3,0,0],[4,0,0],[4,0,0]]
1959                    CSI[1][9:] = [[1./dT,0.,0.],[-dT,0.,0.],[-dT,0.,0.]]
1960                    USC = [1,1,1,1,1,1,1,1,1,1,1,1]                             
1961                elif '(x)' in siteSym:
1962                    CSI[1][9:] = [-dT,0.,0.],[-dT,0.,0.],[1./dT,0.,0.]
1963                elif '(y)' in siteSym:
1964                    CSI[1][9:] = [-dT,0.,0.],[1./dT,0.,0.],[-dT,0.,0.]
1965                elif '(z)' in siteSym:
1966                    CSI[1][9:] = [1./dT,0.,0.],[-dT,0.,0.],[-dT,0.,0.]
1967                for i in range(6):
1968                    if not USC[i]:
1969                        CSI[0][i] = [0,0,0]
1970                        CSI[1][i] = [0.,0.,0.]
1971                        CSI[0][i+6] = [0,0,0]
1972                        CSI[1][i+6] = [0.,0.,0.]
1973            else:                       
1974                for i in range(6):
1975                    if not np.allclose(dU[i,i,:],dUT[i,i,:]*sdet):  #sin part
1976                        usc[i] = 0
1977                    if not np.allclose(dU[i,i+6,:],dUT[i,i+6,:]):   #cos part
1978                        usc[i+6] = 0
1979                if np.any(dUT[1,0,:]):
1980                    if '4/m' in siteSym:
1981                        CSI[0][6:8] = [[12,0,0],[12,0,0]]
1982                        if ssop[1][3]:
1983                            CSI[1][6:8] = [[1.,0.,0.],[-1.,0.,0.]]
1984                            usc[9] = 1
1985                        else:
1986                            CSI[1][6:8] = [[1.,0.,0.],[1.,0.,0.]]
1987                            usc[9] = 0
1988                    elif '4' in siteSym:
1989                        CSI[0][6:8] = [[12,0,0],[12,0,0]]
1990                        CSI[0][:2] = [[11,0,0],[11,0,0]]
1991                        if ssop[1][3]:
1992                            CSI[1][:2] = [[1.,0.,0.],[-1.,0.,0.]]
1993                            CSI[1][6:8] = [[1.,0.,0.],[-1.,0.,0.]]
1994                            usc[2] = 0
1995                            usc[8] = 0
1996                            usc[3] = 1
1997                            usc[9] = 1
1998                        else:
1999                            CSI[1][:2] = [[1.,0.,0.],[1.,0.,0.]]
2000                            CSI[1][6:8] = [[1.,0.,0.],[1.,0.,0.]]
2001                            usc[2] = 1
2002                            usc[8] = 1
2003                            usc[3] = 0               
2004                            usc[9] = 0
2005                    elif 'xy' in siteSym or '+-0' in siteSym:
2006                        if np.allclose(dU[0,0,:],dUT[0,1,:]*sdet):
2007                            CSI[0][4:6] = [[12,0,0],[12,0,0]]
2008                            CSI[0][6:8] = [[11,0,0],[11,0,0]]
2009                            CSI[1][4:6] = [[1.,0.,0.],[sdet,0.,0.]]
2010                            CSI[1][6:8] = [[1.,0.,0.],[sdet,0.,0.]]
2011                            usc[4:6] = 0
2012                            usc[6:8] = 0
2013                       
2014                print SSMT2text(ssop).replace(' ',''),sdet,ssdet,epsinv,usc
2015            USC &= usc
2016        print USC
2017        if not np.any(dtau%.5):
2018            n = -1
2019            for i,U in enumerate(USC):
2020                if U:
2021                    n += 1
2022                    CSI[0][i][0] = n+1
2023                    CSI[1][i][0] = 1.0
2024
2025        return CSI,dU,dUTP
2026       
2027    def DoUiso():
2028        print 'Uiso'
2029        delt4 = np.eye(4)*0.001
2030       
2031    print 'super space group: ',SSGData['SSpGrp']
2032    CSI = {'Sfrac':[[[1,0],[2,0]],[[1.,0.],[1.,0.]]],
2033        'Spos':[[[1,0,0],[2,0,0],[3,0,0], [4,0,0],[5,0,0],[6,0,0]],
2034            [[1.,0.,0.],[1.,0.,0.],[1.,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.]]],    #sin & cos
2035        'Sadp':[[[1,0,0],[2,0,0],[3,0,0],[4,0,0],[5,0,0],[6,0,0], 
2036            [7,0,0],[8,0,0],[9,0,0],[10,0,0],[11,0,0],[12,0,0]],
2037            [[1.,0.,0.],[1.,0.,0.],[1.,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.],
2038            [1.,0.,0.],[1.,0.,0.],[1.,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.]]],
2039        'Smag':[[[1,0,0],[2,0,0],[3,0,0], [4,0,0],[5,0,0],[6,0,0]],
2040            [[1.,0.,0.],[1.,0.,0.],[1.,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.]]],}
2041    xyz = np.array(XYZ)%1.
2042    xyzt = np.array(XYZ+[0,])%1.
2043    SGOps = copy.deepcopy(SGData['SGOps'])
2044    siteSym = SytSym(XYZ,SGData)[0].strip()
2045    print 'siteSym: ',siteSym
2046    if siteSym == '1':   #"1" site symmetry
2047        if debug:
2048            return CSI,None,None,None,None
2049        else:
2050            return CSI
2051    elif siteSym == '-1':   #"-1" site symmetry
2052        CSI['Sfrac'][0] = [[1,0],[0,0]]
2053        CSI['Spos'][0] = [[1,0,0],[2,0,0],[3,0,0], [0,0,0],[0,0,0],[0,0,0]]
2054        CSI['Sadp'][0] = [[0,0,0],[0,0,0],[0,0,0],[0,0,0],[0,0,0],[0,0,0], 
2055        [1,0,0],[2,0,0],[3,0,0],[4,0,0],[5,0,0],[6,0,0]]
2056        if debug:
2057            return CSI,None,None,None,None
2058        else:
2059            return CSI
2060    SSGOps = copy.deepcopy(SSGData['SSGOps'])
2061    #expand ops to include inversions if any
2062    if SGData['SGInv']:
2063        for op,sop in zip(SGData['SGOps'],SSGData['SSGOps']):
2064            SGOps.append([-op[0],-op[1]%1.])
2065            SSGOps.append([-sop[0],-sop[1]%1.])
2066    #build set of sym ops around special position       
2067    SSop = []
2068    Sop = []
2069    Sdtau = []
2070    for iop,Op in enumerate(SGOps):         
2071        nxyz = (np.inner(Op[0],xyz)+Op[1])%1.
2072        if np.allclose(xyz,nxyz,1.e-4) and iop and MT2text(Op).replace(' ','') != '-X,-Y,-Z':
2073            SSop.append(SSGOps[iop])
2074            Sop.append(SGOps[iop])
2075            ssopinv = nl.inv(SSGOps[iop][0])
2076            mst = ssopinv[3][:3]
2077            epsinv = ssopinv[3][3]
2078            Sdtau.append(np.sum(mst*(XYZ-SGOps[iop][1])-epsinv*SSGOps[iop][1][3]))
2079    SdIndx = np.argsort(np.array(Sdtau))     # just to do in sensible order
2080    OpText =  [MT2text(s).replace(' ','') for s in Sop]         #debug?
2081    SSOpText = [SSMT2text(ss).replace(' ','') for ss in SSop]   #debug?
2082    print 'special pos super operators: ',SSOpText
2083    #setup displacement arrays
2084    tau = np.linspace(-1,1,49,True)
2085    #make modulation arrays - one parameter at a time
2086    #site fractions
2087    CSI['Sfrac'],dF,dFTP = DoFrac()
2088    #positions
2089    CSI['Spos'],dX,dXTP = DoXYZ()       
2090    #anisotropic thermal motion
2091    CSI['Sadp'],dU,dUTP = DoUij()
2092    CSI['Spos'][0] = orderParms(CSI['Spos'][0])
2093    CSI['Sadp'][0] = orderParms(CSI['Sadp'][0])           
2094    if debug:
2095        return CSI,tau,[dF,dFTP],[dX,dXTP],[dU,dUTP]
2096    else:
2097        return CSI
2098   
2099def MustrainNames(SGData):
2100    'Needs a doc string'
2101    laue = SGData['SGLaue']
2102    uniq = SGData['SGUniq']
2103    if laue in ['m3','m3m']:
2104        return ['S400','S220']
2105    elif laue in ['6/m','6/mmm','3m1']:
2106        return ['S400','S004','S202']
2107    elif laue in ['31m','3']:
2108        return ['S400','S004','S202','S211']
2109    elif laue in ['3R','3mR']:
2110        return ['S400','S220','S310','S211']
2111    elif laue in ['4/m','4/mmm']:
2112        return ['S400','S004','S220','S022']
2113    elif laue in ['mmm']:
2114        return ['S400','S040','S004','S220','S202','S022']
2115    elif laue in ['2/m']:
2116        SHKL = ['S400','S040','S004','S220','S202','S022']
2117        if uniq == 'a':
2118            SHKL += ['S013','S031','S211']
2119        elif uniq == 'b':
2120            SHKL += ['S301','S103','S121']
2121        elif uniq == 'c':
2122            SHKL += ['S130','S310','S112']
2123        return SHKL
2124    else:
2125        SHKL = ['S400','S040','S004','S220','S202','S022']
2126        SHKL += ['S310','S103','S031','S130','S301','S013']
2127        SHKL += ['S211','S121','S112']
2128        return SHKL
2129       
2130def HStrainVals(HSvals,SGData):
2131    laue = SGData['SGLaue']
2132    uniq = SGData['SGUniq']
2133    DIJ = np.zeros(6)
2134    if laue in ['m3','m3m']:
2135        DIJ[:3] = [HSvals[0],HSvals[0],HSvals[0]]
2136    elif laue in ['6/m','6/mmm','3m1','31m','3']:
2137        DIJ[:4] = [HSvals[0],HSvals[0],HSvals[1],HSvals[0]]
2138    elif laue in ['3R','3mR']:
2139        DIJ = [HSvals[0],HSvals[0],HSvals[0],HSvals[1],HSvals[1],HSvals[1]]
2140    elif laue in ['4/m','4/mmm']:
2141        DIJ[:3] = [HSvals[0],HSvals[0],HSvals[1]]
2142    elif laue in ['mmm']:
2143        DIJ[:3] = [HSvals[0],HSvals[1],HSvals[2]]
2144    elif laue in ['2/m']:
2145        DIJ[:3] = [HSvals[0],HSvals[1],HSvals[2]]
2146        if uniq == 'a':
2147            DIJ[5] = HSvals[3]
2148        elif uniq == 'b':
2149            DIJ[4] = HSvals[3]
2150        elif uniq == 'c':
2151            DIJ[3] = HSvals[3]
2152    else:
2153        DIJ = [HSvals[0],HSvals[1],HSvals[2],HSvals[3],HSvals[4],HSvals[5]]
2154    return DIJ
2155
2156def HStrainNames(SGData):
2157    'Needs a doc string'
2158    laue = SGData['SGLaue']
2159    uniq = SGData['SGUniq']
2160    if laue in ['m3','m3m']:
2161        return ['D11','eA']         #add cubic strain term
2162    elif laue in ['6/m','6/mmm','3m1','31m','3']:
2163        return ['D11','D33']
2164    elif laue in ['3R','3mR']:
2165        return ['D11','D12']
2166    elif laue in ['4/m','4/mmm']:
2167        return ['D11','D33']
2168    elif laue in ['mmm']:
2169        return ['D11','D22','D33']
2170    elif laue in ['2/m']:
2171        Dij = ['D11','D22','D33']
2172        if uniq == 'a':
2173            Dij += ['D23']
2174        elif uniq == 'b':
2175            Dij += ['D13']
2176        elif uniq == 'c':
2177            Dij += ['D12']
2178        return Dij
2179    else:
2180        Dij = ['D11','D22','D33','D12','D13','D23']
2181        return Dij
2182   
2183def MustrainCoeff(HKL,SGData):
2184    'Needs a doc string'
2185    #NB: order of terms is the same as returned by MustrainNames
2186    laue = SGData['SGLaue']
2187    uniq = SGData['SGUniq']
2188    h,k,l = HKL
2189    Strm = []
2190    if laue in ['m3','m3m']:
2191        Strm.append(h**4+k**4+l**4)
2192        Strm.append(3.0*((h*k)**2+(h*l)**2+(k*l)**2))
2193    elif laue in ['6/m','6/mmm','3m1']:
2194        Strm.append(h**4+k**4+2.0*k*h**3+2.0*h*k**3+3.0*(h*k)**2)
2195        Strm.append(l**4)
2196        Strm.append(3.0*((h*l)**2+(k*l)**2+h*k*l**2))
2197    elif laue in ['31m','3']:
2198        Strm.append(h**4+k**4+2.0*k*h**3+2.0*h*k**3+3.0*(h*k)**2)
2199        Strm.append(l**4)
2200        Strm.append(3.0*((h*l)**2+(k*l)**2+h*k*l**2))
2201        Strm.append(4.0*h*k*l*(h+k))
2202    elif laue in ['3R','3mR']:
2203        Strm.append(h**4+k**4+l**4)
2204        Strm.append(3.0*((h*k)**2+(h*l)**2+(k*l)**2))
2205        Strm.append(2.0*(h*l**3+l*k**3+k*h**3)+2.0*(l*h**3+k*l**3+l*k**3))
2206        Strm.append(4.0*(k*l*h**2+h*l*k**2+h*k*l**2))
2207    elif laue in ['4/m','4/mmm']:
2208        Strm.append(h**4+k**4)
2209        Strm.append(l**4)
2210        Strm.append(3.0*(h*k)**2)
2211        Strm.append(3.0*((h*l)**2+(k*l)**2))
2212    elif laue in ['mmm']:
2213        Strm.append(h**4)
2214        Strm.append(k**4)
2215        Strm.append(l**4)
2216        Strm.append(3.0*(h*k)**2)
2217        Strm.append(3.0*(h*l)**2)
2218        Strm.append(3.0*(k*l)**2)
2219    elif laue in ['2/m']:
2220        Strm.append(h**4)
2221        Strm.append(k**4)
2222        Strm.append(l**4)
2223        Strm.append(3.0*(h*k)**2)
2224        Strm.append(3.0*(h*l)**2)
2225        Strm.append(3.0*(k*l)**2)
2226        if uniq == 'a':
2227            Strm.append(2.0*k*l**3)
2228            Strm.append(2.0*l*k**3)
2229            Strm.append(4.0*k*l*h**2)
2230        elif uniq == 'b':
2231            Strm.append(2.0*l*h**3)
2232            Strm.append(2.0*h*l**3)
2233            Strm.append(4.0*h*l*k**2)
2234        elif uniq == 'c':
2235            Strm.append(2.0*h*k**3)
2236            Strm.append(2.0*k*h**3)
2237            Strm.append(4.0*h*k*l**2)
2238    else:
2239        Strm.append(h**4)
2240        Strm.append(k**4)
2241        Strm.append(l**4)
2242        Strm.append(3.0*(h*k)**2)
2243        Strm.append(3.0*(h*l)**2)
2244        Strm.append(3.0*(k*l)**2)
2245        Strm.append(2.0*k*h**3)
2246        Strm.append(2.0*h*l**3)
2247        Strm.append(2.0*l*k**3)
2248        Strm.append(2.0*h*k**3)
2249        Strm.append(2.0*l*h**3)
2250        Strm.append(2.0*k*l**3)
2251        Strm.append(4.0*k*l*h**2)
2252        Strm.append(4.0*h*l*k**2)
2253        Strm.append(4.0*k*h*l**2)
2254    return Strm
2255   
2256def Muiso2Shkl(muiso,SGData,cell):
2257    "this is to convert isotropic mustrain to generalized Shkls"
2258    import GSASIIlattice as G2lat
2259    A = G2lat.cell2AB(cell)[0]
2260   
2261    def minMus(Shkl,muiso,H,SGData,A):
2262        U = np.inner(A.T,H)
2263        S = np.array(MustrainCoeff(U,SGData))
2264        Sum = np.sqrt(np.sum(np.multiply(S,Shkl[:,np.newaxis]),axis=0))
2265        rad = np.sqrt(np.sum((Sum[:,np.newaxis]*H)**2,axis=1))
2266        return (muiso-rad)**2
2267       
2268    laue = SGData['SGLaue']
2269    PHI = np.linspace(0.,360.,60,True)
2270    PSI = np.linspace(0.,180.,60,True)
2271    X = np.outer(npsind(PHI),npsind(PSI))
2272    Y = np.outer(npcosd(PHI),npsind(PSI))
2273    Z = np.outer(np.ones(np.size(PHI)),npcosd(PSI))
2274    HKL = np.dstack((X,Y,Z))
2275    if laue in ['m3','m3m']:
2276        S0 = [1000.,1000.]
2277    elif laue in ['6/m','6/mmm','3m1']:
2278        S0 = [1000.,1000.,1000.]
2279    elif laue in ['31m','3']:
2280        S0 = [1000.,1000.,1000.,1000.]
2281    elif laue in ['3R','3mR']:
2282        S0 = [1000.,1000.,1000.,1000.]
2283    elif laue in ['4/m','4/mmm']:
2284        S0 = [1000.,1000.,1000.,1000.]
2285    elif laue in ['mmm']:
2286        S0 = [1000.,1000.,1000.,1000.,1000.,1000.]
2287    elif laue in ['2/m']:
2288        S0 = [1000.,1000.,1000.,0.,0.,0.,0.,0.,0.]
2289    else:
2290        S0 = [1000.,1000.,1000.,1000.,1000., 1000.,1000.,1000.,1000.,1000., 
2291            1000.,1000.,0.,0.,0.]
2292    S0 = np.array(S0)
2293    HKL = np.reshape(HKL,(-1,3))
2294    result = so.leastsq(minMus,S0,(np.ones(HKL.shape[0])*muiso,HKL,SGData,A))
2295    return result[0]
2296       
2297def SytSym(XYZ,SGData):
2298    '''
2299    Generates the number of equivalent positions and a site symmetry code for a specified coordinate and space group
2300
2301    :param XYZ: an array, tuple or list containing 3 elements: x, y & z
2302    :param SGData: from SpcGroup
2303    :Returns: a two element tuple:
2304
2305     * The 1st element is a code for the site symmetry (see GetKNsym)
2306     * The 2nd element is the site multiplicity
2307
2308    '''
2309    def PackRot(SGOps):
2310        IRT = []
2311        for ops in SGOps:
2312            M = ops[0]
2313            irt = 0
2314            for j in range(2,-1,-1):
2315                for k in range(2,-1,-1):
2316                    irt *= 3
2317                    irt += M[k][j]
2318            IRT.append(int(irt))
2319        return IRT
2320       
2321    SymName = ''
2322    Mult = 1
2323    Isym = 0
2324    if SGData['SGLaue'] in ['3','3m1','31m','6/m','6/mmm']:
2325        Isym = 1073741824
2326    Jdup = 0
2327    Xeqv = GenAtom(XYZ,SGData,True)
2328    IRT = PackRot(SGData['SGOps'])
2329    L = -1
2330    for ic,cen in enumerate(SGData['SGCen']):
2331        for invers in range(int(SGData['SGInv']+1)):
2332            for io,ops in enumerate(SGData['SGOps']):
2333                irtx = (1-2*invers)*IRT[io]
2334                L += 1
2335                if not Xeqv[L][1]:
2336                    Jdup += 1
2337                    jx = GetOprPtrName(str(irtx))
2338                    if jx[2] < 39:
2339                        Isym += 2**(jx[2]-1)
2340    if Isym == 1073741824: Isym = 0
2341    Mult = len(SGData['SGOps'])*len(SGData['SGCen'])*(int(SGData['SGInv'])+1)/Jdup
2342         
2343    return GetKNsym(str(Isym)),Mult
2344   
2345def ElemPosition(SGData):
2346    ''' Under development.
2347    Object here is to return a list of symmetry element types and locations suitable
2348    for say drawing them.
2349    So far I have the element type... getting all possible locations without lookup may be impossible!
2350    '''
2351    SymElements = []
2352    Inv = SGData['SGInv']
2353    Cen = SGData['SGCen']
2354    eleSym = {-3:['','-1'],-2:['',-6],-1:['2','-4'],0:['3','-3'],1:['4','m'],2:['6',''],3:['1','']}
2355    # get operators & expand if centrosymmetric
2356    Ops = SGData['SGOps']
2357    opM = np.array([op[0].T for op in Ops])
2358    opT = np.array([op[1] for op in Ops])
2359    if Inv:
2360        opM = np.concatenate((opM,-opM))
2361        opT = np.concatenate((opT,-opT))
2362    opMT = zip(opM,opT)
2363    for M,T in opMT[1:]:        #skip I
2364        Dt = int(nl.det(M))
2365        Tr = int(np.trace(M))
2366        Dt = -(Dt-1)/2
2367        Es = eleSym[Tr][Dt]
2368        if Dt:              #rotation-inversion
2369            I = np.eye(3)
2370            if Tr == 1:     #mirrors/glides
2371                if np.any(T):       #glide
2372                    M2 = np.inner(M,M)
2373                    MT = np.inner(M,T)+T
2374                    print 'glide',Es,MT
2375                    print M2
2376                else:               #mirror
2377                    print 'mirror',Es,T
2378                    print I-M
2379                X = [-1,-1,-1]
2380            elif Tr == -3:  # pure inversion
2381                X = np.inner(nl.inv(I-M),T)
2382                print 'inversion',Es,X
2383            else:           #other rotation-inversion
2384                M2 = np.inner(M,M)
2385                MT = np.inner(M,T)+T
2386                print 'rot-inv',Es,MT
2387                print M2
2388                X = [-1,-1,-1]
2389        else:               #rotations
2390            print 'rotation',Es
2391            X = [-1,-1,-1]
2392        #SymElements.append([Es,X])
2393       
2394    return #SymElements
2395   
2396def ApplyStringOps(A,SGData,X,Uij=[]):
2397    'Needs a doc string'
2398    SGOps = SGData['SGOps']
2399    SGCen = SGData['SGCen']
2400    Ax = A.split('+')
2401    Ax[0] = int(Ax[0])
2402    iC = 0
2403    if Ax[0] < 0:
2404        iC = 1
2405    Ax[0] = abs(Ax[0])
2406    nA = Ax[0]%100-1
2407    cA = Ax[0]/100
2408    Cen = SGCen[cA]
2409    M,T = SGOps[nA]
2410    if len(Ax)>1:
2411        cellA = Ax[1].split(',')
2412        cellA = np.array([int(a) for a in cellA])
2413    else:
2414        cellA = np.zeros(3)
2415    newX = Cen+(1-2*iC)*(np.inner(M,X)+T)+cellA
2416    if len(Uij):
2417        U = Uij2U(Uij)
2418        U = np.inner(M,np.inner(U,M).T)
2419        newUij = U2Uij(U)
2420        return [newX,newUij]
2421    else:
2422        return newX
2423       
2424def StringOpsProd(A,B,SGData):
2425    """
2426    Find A*B where A & B are in strings '-' + '100*c+n' + '+ijk'
2427    where '-' indicates inversion, c(>0) is the cell centering operator,
2428    n is operator number from SgOps and ijk are unit cell translations (each may be <0).
2429    Should return resultant string - C. SGData - dictionary using entries:
2430
2431       *  'SGCen': cell centering vectors [0,0,0] at least
2432       *  'SGOps': symmetry operations as [M,T] so that M*x+T = x'
2433
2434    """
2435    SGOps = SGData['SGOps']
2436    SGCen = SGData['SGCen']
2437    #1st split out the cell translation part & work on the operator parts
2438    Ax = A.split('+'); Bx = B.split('+')
2439    Ax[0] = int(Ax[0]); Bx[0] = int(Bx[0])
2440    iC = 0
2441    if Ax[0]*Bx[0] < 0:
2442        iC = 1
2443    Ax[0] = abs(Ax[0]); Bx[0] = abs(Bx[0])
2444    nA = Ax[0]%100-1;  nB = Bx[0]%100-1
2445    cA = Ax[0]/100;  cB = Bx[0]/100
2446    Cen = (SGCen[cA]+SGCen[cB])%1.0
2447    cC = np.nonzero([np.allclose(C,Cen) for C in SGCen])[0][0]
2448    Ma,Ta = SGOps[nA]; Mb,Tb = SGOps[nB]
2449    Mc = np.inner(Ma,Mb.T)
2450#    print Ma,Mb,Mc
2451    Tc = (np.add(np.inner(Mb,Ta)+1.,Tb))%1.0
2452#    print Ta,Tb,Tc
2453#    print [np.allclose(M,Mc)&np.allclose(T,Tc) for M,T in SGOps]
2454    nC = np.nonzero([np.allclose(M,Mc)&np.allclose(T,Tc) for M,T in SGOps])[0][0]
2455    #now the cell translation part
2456    if len(Ax)>1:
2457        cellA = Ax[1].split(',')
2458        cellA = [int(a) for a in cellA]
2459    else:
2460        cellA = [0,0,0]
2461    if len(Bx)>1:
2462        cellB = Bx[1].split(',')
2463        cellB = [int(b) for b in cellB]
2464    else:
2465        cellB = [0,0,0]
2466    cellC = np.add(cellA,cellB)
2467    C = str(((nC+1)+(100*cC))*(1-2*iC))+'+'+ \
2468        str(int(cellC[0]))+','+str(int(cellC[1]))+','+str(int(cellC[2]))
2469    return C
2470           
2471def U2Uij(U):
2472    #returns the UIJ vector U11,U22,U33,U12,U13,U23 from tensor U
2473    return [U[0][0],U[1][1],U[2][2],2.*U[0][1],2.*U[0][2],2.*U[1][2]]
2474   
2475def Uij2U(Uij):
2476    #returns the thermal motion tensor U from Uij as numpy array
2477    return np.array([[Uij[0],Uij[3]/2.,Uij[4]/2.],[Uij[3]/2.,Uij[1],Uij[5]/2.],[Uij[4]/2.,Uij[5]/2.,Uij[2]]])
2478
2479def StandardizeSpcName(spcgroup):
2480    '''Accept a spacegroup name where spaces may have not been used
2481    in the names according to the GSAS convention (spaces between symmetry
2482    for each axis) and return the space group name as used in GSAS
2483    '''
2484    rspc = spcgroup.replace(' ','').upper()
2485    # deal with rhombohedral and hexagonal setting designations
2486    rhomb = ''
2487    if rspc[-1:] == 'R':
2488        rspc = rspc[:-1]
2489        rhomb = ' R'
2490    elif rspc[-1:] == 'H': # hexagonal is assumed and thus can be ignored
2491        rspc = rspc[:-1]
2492    # look for a match in the spacegroup lists
2493    for i in spglist.values():
2494        for spc in i:
2495            if rspc == spc.replace(' ','').upper():
2496                return spc + rhomb
2497    # how about the post-2002 orthorhombic names?
2498    for i,spc in sgequiv_2002_orthorhombic:
2499        if rspc == i.replace(' ','').upper():
2500            return spc
2501    # not found
2502    return ''
2503
2504   
2505spglist = {}
2506'''A dictionary of space groups as ordered and named in the pre-2002 International
2507Tables Volume A, except that spaces are used following the GSAS convention to
2508separate the different crystallographic directions.
2509Note that the symmetry codes here will recognize many non-standard space group
2510symbols with different settings. They are ordered by Laue group
2511'''
2512spglist = {
2513    'P1' : ('P 1','P -1',), # 1-2
2514    'P2/m': ('P 2','P 21','P m','P a','P c','P n',
2515        'P 2/m','P 21/m','P 2/c','P 2/a','P 2/n','P 21/c','P 21/a','P 21/n',), #3-15
2516    'C2/m':('C 2','C m','C c','C n',
2517        'C 2/m','C 2/c','C 2/n',),
2518    'Pmmm':('P 2 2 2',
2519        'P 2 2 21','P 21 2 2','P 2 21 2',
2520        'P 21 21 2','P 2 21 21','P 21 2 21',
2521        'P 21 21 21',
2522        'P m m 2','P 2 m m','P m 2 m',
2523        'P m c 21','P 21 m a','P b 21 m','P m 21 b','P c m 21','P 21 a m',
2524        'P c c 2','P 2 a a','P b 2 b',
2525        'P m a 2','P 2 m b','P c 2 m','P m 2 a','P b m 2','P 2 c m',
2526        'P c a 21','P 21 a b','P c 21 b','P b 21 a','P b c 21','P 21 c a',
2527        'P n c 2','P 2 n a','P b 2 n','P n 2 b','P c n 2','P 2 a n',
2528        'P m n 21','P 21 m n','P n 21 m','P m 21 n','P n m 21','P 21 n m',
2529        'P b a 2','P 2 c b','P c 2 a',
2530        'P n a 21','P 21 n b','P c 21 n','P n 21 a','P b n 21','P 21 c n',
2531        'P n n 2','P 2 n n','P n 2 n',
2532        'P m m m','P n n n',
2533        'P c c m','P m a a','P b m b',
2534        'P b a n','P n c b','P c n a',
2535        'P m m a','P b m m','P m c m','P m a m','P m m b','P c m m',
2536        'P n n a','P b n n','P n c n','P n a n','P n n b','P c n n',
2537        'P m n a','P b m n','P n c m','P m a n','P n m b','P c n m',
2538        'P c c a','P b a a','P b c b','P b a b','P c c b','P c a a',
2539        'P b a m','P m c b','P c m a',
2540        'P c c n','P n a a','P b n b',
2541        'P b c m','P m c a','P b m a','P c m b','P c a m','P m a b',
2542        'P n n m','P m n n','P n m n',
2543        'P m m n','P n m m','P m n m',
2544        'P b c n','P n c a','P b n a','P c n b','P c a n','P n a b',
2545        'P b c a','P c a b',
2546        'P n m a','P b n m','P m c n','P n a m','P m n b','P c m n',
2547        ),
2548    'Cmmm':('C 2 2 21','C 2 2 2','C m m 2',
2549        'C m c 21','C c m 21','C c c 2','C m 2 m','C 2 m m',
2550        'C m 2 a','C 2 m b','C c 2 m','C 2 c m','C c 2 a','C 2 c b',
2551        'C m c m','C m c a','C c m b',
2552        'C m m m','C c c m','C m m a','C m m b','C c c a','C c c b',),
2553    'Immm':('I 2 2 2','I 21 21 21',
2554        'I m m 2','I m 2 m','I 2 m m',
2555        'I b a 2','I 2 c b','I c 2 a',
2556        'I m a 2','I 2 m b','I c 2 m','I m 2 a','I b m 2','I 2 c m',
2557        'I m m m','I b a m','I m c b','I c m a',
2558        'I b c a','I c a b',
2559        'I m m a','I b m m ','I m c m','I m a m','I m m b','I c m m',),
2560    'Fmmm':('F 2 2 2','F m m m', 'F d d d',
2561        'F m m 2','F m 2 m','F 2 m m',
2562        'F d d 2','F d 2 d','F 2 d d',),
2563    'P4/mmm':('P 4','P 41','P 42','P 43','P -4','P 4/m','P 42/m','P 4/n','P 42/n',
2564        'P 4 2 2','P 4 21 2','P 41 2 2','P 41 21 2','P 42 2 2',
2565        'P 42 21 2','P 43 2 2','P 43 21 2','P 4 m m','P 4 b m','P 42 c m',
2566        'P 42 n m','P 4 c c','P 4 n c','P 42 m c','P 42 b c','P -4 2 m',
2567        'P -4 2 c','P -4 21 m','P -4 21 c','P -4 m 2','P -4 c 2','P -4 b 2',
2568        'P -4 n 2','P 4/m m m','P 4/m c c','P 4/n b m','P 4/n n c','P 4/m b m',
2569        'P 4/m n c','P 4/n m m','P 4/n c c','P 42/m m c','P 42/m c m',
2570        'P 42/n b c','P 42/n n m','P 42/m b c','P 42/m n m','P 42/n m c',
2571        'P 42/n c m',),
2572    'I4/mmm':('I 4','I 41','I -4','I 4/m','I 41/a','I 4 2 2','I 41 2 2','I 4 m m',
2573        'I 4 c m','I 41 m d','I 41 c d',
2574        'I -4 m 2','I -4 c 2','I -4 2 m','I -4 2 d','I 4/m m m','I 4/m c m',
2575        'I 41/a m d','I 41/a c d'),
2576    'R3-H':('R 3','R -3','R 3 2','R 3 m','R 3 c','R -3 m','R -3 c',),
2577    'P6/mmm': ('P 3','P 31','P 32','P -3','P 3 1 2','P 3 2 1','P 31 1 2',
2578        'P 31 2 1','P 32 1 2','P 32 2 1', 'P 3 m 1','P 3 1 m','P 3 c 1',
2579        'P 3 1 c','P -3 1 m','P -3 1 c','P -3 m 1','P -3 c 1','P 6','P 61',
2580        'P 65','P 62','P 64','P 63','P -6','P 6/m','P 63/m','P 6 2 2',
2581        'P 61 2 2','P 65 2 2','P 62 2 2','P 64 2 2','P 63 2 2','P 6 m m',
2582        'P 6 c c','P 63 c m','P 63 m c','P -6 m 2','P -6 c 2','P -6 2 m',
2583        'P -6 2 c','P 6/m m m','P 6/m c c','P 63/m c m','P 63/m m c',),
2584    'Pm3m': ('P 2 3','P 21 3','P m 3','P n 3','P a 3','P 4 3 2','P 42 3 2',
2585        'P 43 3 2','P 41 3 2','P -4 3 m','P -4 3 n','P m 3 m','P n 3 n',
2586        'P m 3 n','P n 3 m',),
2587    'Im3m':('I 2 3','I 21 3','I m -3','I a -3', 'I 4 3 2','I 41 3 2',
2588        'I -4 3 m', 'I -4 3 d','I m -3 m','I m 3 m','I a -3 d',),
2589    'Fm3m':('F 2 3','F m -3','F d -3','F 4 3 2','F 41 3 2','F -4 3 m',
2590        'F -4 3 c','F m -3 m','F m 3 m','F m -3 c','F d -3 m','F d -3 c',),
2591}
2592
2593ssdict = {}
2594'''A dictionary of superspace group symbols allowed for each entry in spglist
2595(except cubics). Monoclinics are all b-unique setting.
2596'''
2597ssdict = {
2598#1,2
2599    'P 1':['(abg)',],'P -1':['(abg)',],
2600#monoclinic - done
2601#3
2602    'P 2':['(a0g)','(a1/2g)','(0b0)','(0b0)s','(1/2b0)','(0b1/2)',],
2603#4       
2604    'P 21':['(a0g)','(0b0)','(1/2b0)','(0b1/2)',],
2605#5
2606    'C 2':['(a0g)','(0b0)','(0b0)s','(0b1/2)',],
2607#6
2608    'P m':['(a0g)','(a0g)s','(a1/2g)','(0b0)','(1/2b0)','(0b1/2)',],
2609#7
2610    'P a':['(a0g)','(a1/2g)','(0b0)','(0b1/2)',],
2611    'P c':['(a0g)','(a1/2g)','(0b0)','(1/2b0)',],
2612    'P n':['(a0g)','(a1/2g)','(0b0)','(1/2b1/2)',],
2613#8       
2614    'C m':['(a0g)','(a0g)s','(0b0)','(0b1/2)',],
2615#9       
2616    'C c':['(a0g)','(a0g)s','(0b0)',],
2617    'C n':['(a0g)','(a0g)s','(0b0)',],
2618#10       
2619    'P 2/m':['(a0g)','(a0g)0s','(a1/2g)','(0b0)','(0b0)s0','(1/2b0)','(0b1/2)',],
2620#11
2621    'P 21/m':['(a0g)','(a0g)0s','(0b0)','(0b0)s0','(1/2b0)','(0b1/2)',],
2622#12       
2623    'C 2/m':['(a0g)','(a0g)0s','(0b0)','(0b0)s0','(0b1/2)',],
2624#13
2625    'P 2/c':['(a0g)','(a0g)0s','(a1/2g)','(0b0)','(0b0)s0','(1/2b0)',],
2626    'P 2/a':['(a0g)','(a0g)0s','(a1/2g)','(0b0)','(0b0)s0','(0b1/2)',],
2627    'P 2/n':['(a0g)','(a0g)0s','(a1/2g)','(0b0)','(0b0)s0','(1/2b1/2)',],
2628#14
2629    'P 21/c':['(a0g)','(0b0)','(1/2b0)',],
2630    'P 21/a':['(a0g)','(0b0)','(0b1/2)',],
2631    'P 21/n':['(a0g)','(0b0)','(1/2b1/2)',],
2632#15
2633    'C 2/c':['(a0g)','(0b0)','(0b0)s0',],
2634    'C 2/n':['(a0g)','(0b0)','(0b0)s0',],
2635#orthorhombic
2636#16   
2637    'P 2 2 2':['(00g)','(00g)00s','(01/2g)','(1/20g)','(1/21/2g)',
2638        '(a00)','(a00)s00','(a01/2)','(a1/20)','(a1/21/2)',
2639        '(0b0)','(0b0)0s0','(1/2b0)','(0b1/2)','(1/2b1/2)',],
2640#17       
2641    'P 2 2 21':['(00g)','(01/2g)','(1/20g)','(1/21/2g)',
2642        '(a00)','(a00)s00','(a1/20)','(0b0)','(0b0)0s0','(1/2b0)',],
2643    'P 21 2 2':['(a00)','(a01/2)','(a1/20)','(a1/21/2)',
2644        '(0b0)','(0b0)0s0','(1/2b0)','(00g)','(00g)00s','(1/20g)',],
2645    'P 2 21 2':['(0b0)','(0b1/2)','(1/2b0)','(1/2b1/2)',
2646        '(00g)','(00g)00s','(1/20g)','(a00)','(a00)s00','(a1/20)',],
2647#18       
2648    'P 21 21 2':['(00g)','(00g)00s','(a00)','(a01/2)','(0b0)','(0b1/2)',],
2649    'P 2 21 21':['(a00)','(a00)s00','(0b0)','(0b1/2)','(00g)','(01/2g)',],
2650    'P 21 2 21':['(0b0)','(0b0)0s0','(00g)','(01/2g)','(a00)','(a01/2)',],
2651#19       
2652    'P 21 21 21':['(00g)','(a00)','(0b0)',],
2653#20       
2654    'C 2 2 21':['(00g)','(10g)','(01g)','(a00)','(a00)s00','(0b0)','(0b0)0s0',],
2655    'A 21 2 2':['(a00)','(a10)','(a01)','(0b0)','(0b0)0s0','(00g)','(00g)00s',],
2656    'B 2 21 2':['(0b0)','(1b0)','(0b1)','(00g)','(00g)00s','(a00)','(a00)s00',],
2657#21       
2658    'C 2 2 2':['(00g)','(00g)00s','(10g)','(10g)00s','(01g)','(01g)00s',
2659        '(a00)','(a00)s00','(a01/2)','(0b0)','(0b0)0s0','(0b1/2)',],
2660    'A 2 2 2':['(a00)','(a00)s00','(a10)','(a10)s00','(a01)','(a01)s00',
2661        '(0b0)','(0b0)0s0','(1/2b0)','(00g)','(00g)00s','(1/20g)',],
2662    'B 2 2 2':['(0b0)','(0b0)0s0','(1b0)','(1b0)0s0','(0b1)','(0b1)0s0',
2663        '(00g)','(00g)00s','(01/2g)','(a00)','(a00)s00','(a1/20)',],
2664#22       
2665    'F 2 2 2':['(00g)','(00g)00s','(10g)','(01g)',
2666        '(a00)','(a00)s00','(a10)','(a01)',
2667        '(0b0)','(0b0)0s0','(1b0)','(0b1)',],
2668#23       
2669    'I 2 2 2':['(00g)','(00g)00s','(a00)','(a00)s00','(0b0)','(0b0)0s0',],
2670#24       
2671    'I 21 21 21':['(00g)','(00g)00s','(a00)','(a00)s00','(0b0)','(0b0)0s0',],
2672#25       
2673    'P m m 2':['(00g)','(00g)s0s','(00g)0ss','(00g)ss0',
2674        '(01/2g)','(01/2g)s0s','(1/20g)','(1/20g)0ss','(1/21/2g)',
2675        '(a00)','(a00)0s0','(a1/20)','(a01/2)','(a01/2)0s0','(a1/21/2)',
2676        '(0b0)','(0b0)s00','(0b1/2)','(0b1/2)s00','(1/2b0)','(1/2b1/2)',],       
2677    'P 2 m m':['(a00)','(a00)ss0','(a00)s0s','(a00)0ss',
2678        '(a01/2)','(a01/2)ss0','(a1/20)','(a1/20)s0s','(a1/21/2)',
2679        '(0b0)','(0b0)00s','(1/2b0)','(0b1/2)','(0b1/2)00s','(1/2b1/2)',
2680        '(00g)','(00g)0s0','(01/2g)','(01/2g)0s0','(1/20g)','(1/21/2g)',],
2681    'P m 2 m':['(0b0)','(0b0)ss0','(0b0)0ss','(0b0)s0s',
2682        '(0b1/2)','(0b1/2)ss0','(1/2b0)','(1/2b0)0ss','(1/2b1/2)',
2683        '(00g)','(00g)s00','(1/20g)','(01/2g)','(01/2g)s00','(1/21/2g)',
2684        '(a00)','(a00)0s0','(a01/2)','(a01/2)0s0','(a1/20)','(a1/21/2)',],       
2685#26       
2686    'P m c 21':['(00g)','(00g)s0s','(01/2g)','(01/2g)s0s','(1/20g)','(1/21/2g)',
2687        '(a00)','(a00)0s0','(a1/20)','(0b0)','(0b0)s00','(0b1/2)',],
2688    'P 21 m a':['(a00)','(a00)ss0','(a01/2)','(a01/2)ss0','(a1/20)','(a1/21/2)',
2689        '(0b0)','(0b0)00s','(1/2b0)','(00g)','(00g)0s0','(01/2g)',],
2690    'P b 21 m':['(0b0)','(0b0)ss0','(0b1/2)','(0b1/2)ss0','(1/2b0)','(1/2b1/2)',
2691        '(00g)','(00g)s00','(1/20g)','(a00)','(a00)0s0','(a01/2)',],
2692    'P m 21 b':['(a00)','(a00)ss0','(a01/2)','(a01/2)ss0','(a1/20)','(a1/21/2)',
2693        '(00g)','(00g)0s0','(01/2g)','(0b0)','(0b0)s00','(0b1/2)',],
2694    'P c m 21':['(00g)','(00g)0ss','(1/20g)','(1/20g)0ss','(01/2g)','(1/21/2g)',
2695        '(0b0)','(0b0)s00','(1/2b0)','(a00)','(a00)0s0','(a01/2)',],
2696    'P 21 a m':['(0b0)','(0b0)ss0','(0b1/2)','(0b1/2)ss0','(1/2b0)','(1/2b1/2)',
2697        '(a00)','(a00)00s','(a1/20)','(00g)','(00g)s00','(1/20g)',],
2698#27       
2699    'P c c 2':['(00g)','(00g)s0s','(00g)0ss','(01/2g)','(1/20g)','(1/21/2g)',
2700        '(a00)','(a00)0s0','(a1/20)','(0b0)','(0b0)s00','(1/2b0)',],
2701    'P 2 a a':['(a00)','(a00)ss0','(a00)s0s','(a01/2)','(a1/20)','(a1/21/2)',
2702        '(0b0)','(0b0)00s','(0b1/2)','(00g)','(00g)0s0','(01/2g)',],
2703    'P b 2 b':['(0b0)','(0b0)0ss','(0b0)ss0','(1/2b0)','(0b1/2)','(1/2b1/2)',
2704        '(00g)','(00g)s00','(1/20g)','(a00)','(a00)00s','(a01/2)',],
2705#28       
2706    'P m a 2':['(00g)','(00g)s0s','(00g)ss0','(00g)0ss','(01/2g)','(01/2g)s0s',
2707        '(0b1/2)','(0b1/2)s00','(a01/2)','(a00)','(0b0)','(0b0)0s0','(a1/20)','(a1/21/2)'],
2708    'P 2 m b':['(a00)','(a00)s0s','(a00)ss0','(a00)0ss','(a01/2)','(a01/2)s0s',
2709        '(1/20g)','(1/20g)s00','(1/2b0)','(0b0)','(00g)','(00g)0s0','(0b1/2)','(1/2b1/2)'],
2710    'P c 2 m':['(0b0)','(0b0)s0s','(0b0)ss0','(0b0)0ss','(1/2b0)','(1/2b0)s0s',
2711        '(a1/20)','(a1/20)s00','(01/2g)','(00g)','(a00)','(a00)0s0','(1/20g)','(1/21/2g)'],
2712    'P m 2 a':['(0b0)','(0b0)s0s','(0b0)ss0','(0b0)0ss','(0b1/2)','(0b1/2)s0s',
2713        '(01/2g)','(01/2g)s00','(a1/20)','(a00)','(00g)','(00g)0s0','(a01/2)','(a1/21/2)'],
2714    'P b m 2':['(00g)','(00g)s0s','(00g)ss0','(00g)0ss','(1/20g)','(1/20g)s0s',
2715        '(a01/2)','(a01/2)s00','(0b1/2)','(0b0)','(a00)','(a00)0s0','(1/2b0)','(1/2b1/2)'],
2716    'P 2 c m':['(a00)','(a00)s0s','(a00)ss0','(a00)0ss','(a1/20)','(a1/20)s0s',
2717        '(1/2b0)','(1/2b0)s00','(1/20g)','(00g)','(0b0)','(0b0)0s0','(01/2g)','(1/21/2g)'],
2718#29       
2719    'P c a 21':['(00g)','(00g)0ss','(01/2g)','(1/20g)',
2720        '(a00)','(a00)0s0','(a1/20)','(0b0)','(0b0)s00','(1/2b0)',],
2721    'P 21 a b':['(a00)','(a00)s0s','(a01/2)','(a1/20)',
2722        '(0b0)','(0b0)00s','(0b1/2)','(00g)','(00g)0s0','(01/2g)',],
2723    'P c 21 b':['(0b0)','(0b0)ss0','(1/2b0)','(0b1/2)',
2724        '(00g)','(00g)s00','(1/20g)','(a00)','(a00)00s','(a01/2)',],
2725    'P b 21 a':['(0b0)','(0b0)0ss','(0b1/2)','(1/2b0)',
2726        '(a00)','(a00)00s','(a1/20)','(00g)','(00g)s00','(1/20g)',],
2727    'P b c 21':['(00g)','(00g)s0s','(1/20g)','(01/2g)',
2728        '(0b0)','(0b0)s00','(0b1/2)','(a00)','(a00)0s0','(a1/20)',],
2729    'P 21 c a':['(a00)','(a00)ss0','(a1/20)','(a01/2)',
2730        '(00g)','(00g)0s0','(1/20g)','(0b0)','(0b0)00s','(0b1/2)',],
2731#30       
2732    'P c n 2':['(00g)','(00g)s0s','(01/2g)','(a00)','(0b0)','(0b0)s00',
2733        '(a1/20)','(1/2b1/2)q00',],
2734    'P 2 a n':['(a00)','(a00)ss0','(a01/2)','(0b0)','(00g)','(00g)0s0',
2735        '(0b1/2)','(1/21/2g)0q0',],
2736    'P n 2 b':['(0b0)','(0b0)0ss','(1/2b0)','(00g)','(a00)','(a00)00s',
2737        '(1/20g)','(a1/21/2)00q',],
2738    'P b 2 n':['(0b0)','(0b0)ss0','(0b1/2)','(a00)','(00g)','(00g)s00',
2739        '(a01/2)','(1/21/2g)0ss',],
2740    'P n c 2':['(00g)','(00g)0ss','(1/20g)','(0b0)','(a00)','(a00)0s0',
2741        '(1/2b0)','(a1/21/2)s0s',],
2742    'P 2 n a':['(a00)','(a00)s0s','(a1/20)','(00g)','(0b0)','(0b0)00s',
2743        '(01/2g)','(1/2b1/2)ss0',],
2744#31       
2745    'P m n 21':['(00g)','(00g)s0s','(01/2g)','(01/2g)s0s','(a00)','(0b0)',
2746        '(0b0)s00','(a1/20)',],
2747    'P 21 m n':['(a00)','(a00)ss0','(a01/2)','(a01/2)ss0','(0b0)','(00g)',
2748        '(00g)0s0','(0b1/2)',],
2749    'P n 21 m':['(0b0)','(0b0)0ss','(1/2b0)','(1/2b0)0ss','(00g)','(a00)',
2750        '(a00)00s','(1/20g)',],
2751    'P m 21 n':['(0b0)','(0b0)ss0','(0b1/2)','(0b1/2)ss0','(a00)','(00g)',
2752        '(00g)s00','(a01/2)',],
2753    'P n m 21':['(00g)','(00g)0ss','(1/20g)','(1/20g)0ss','(0b0)','(a00)',
2754        '(a00)0s0','(1/2b0)',],
2755    'P 21 n m':['(a00)','(a00)s0s','(a1/20)','(a1/20)s0s','(00g)','(0b0)',
2756        '(0b0)00s','(01/2g)',],
2757#32       
2758    'P b a 2':['(00g)','(00g)s0s','(00g)0ss','(00g)ss0','(1/21/2g)qq0',
2759        '(a00)','(a01/2)','(0b0)','(0b1/2)',],
2760    'P 2 c b':['(a00)','(a00)ss0','(a00)s0s','(a00)0ss','(a1/21/2)0qq',
2761        '(0b0)','(1/2b0)','(00g)','(1/20g)',],
2762    'P c 2 a':['(0b0)','(0b0)ss0','(0b0)0ss','(0b0)s0s','(1/2b1/2)q0q',
2763        '(00g)','01/2g)','(a00)','(a1/20)',],
2764#33       
2765    'P b n 21':['(00g)','(00g)s0s','(1/21/2g)qq0','(a00)','(0b0)',],
2766    'P 21 c n':['(a00)','(a00)ss0','(a1/21/2)0qq','(0b0)','(00g)',],
2767    'P n 21 a':['(0b0)','(0b0)0ss','(1/2b1/2)q0q','(00g)','(a00)',],
2768    'P c 21 n':['(0b0)','(0b0)ss0','(1/2b1/2)q0q','(a00)','(00g)',],
2769    'P n a 21':['(00g)','(00g)0ss','(1/21/2g)qq0','(0b0)','(a00)',],
2770    'P 21 n b':['(a00)','(a00)s0s','(a1/21/2)0qq','(00g)','(0b0)',],
2771#34       
2772    'P n n 2':['(00g)','(00g)s0s','(00g)0ss','(1/21/2g)qq0',
2773        '(a00)','(a1/21/2)0q0','(a1/21/2)00q','(0b0)','(1/2b1/2)q00','(1/2b1/2)00q',],
2774    'P 2 n n':['(a00)','(a00)ss0','(a00)s0s','(a1/21/2)0qq',
2775        '(0b0)','(1/2b1/2)q00','(1/2b1/2)00q','(00g)','(1/21/2g)0q0','(1/21/2g)q00',],
2776    'P n 2 n':['(0b0)','(0b0)ss0','(0b0)0ss','(1/2b1/2)q0q',
2777        '(00g)','(1/21/2g)0q0','(1/21/2g)q00','(a00)','(a1/21/2)00q','(a1/21/2)0q0',],
2778#35       
2779    'C m m 2':['(00g)','(00g)s0s','(00g)ss0','(10g)','(10g)s0s','(10g)ss0',
2780        '(0b0)','(0b0)s00','(0b1/2)','(0b1/2)s00',],
2781    'A 2 m m':['(a00)','(a00)ss0','(a00)0ss','(a10)','(a10)ss0','(a10)0ss',
2782        '(00g)','(00g)0s0','(1/20g)','(1/20g)0s0',],
2783    'B m 2 m':['(0b0)','(0b0)0ss','(0b0)s0s','(0b1)','(0b1)0ss','(0b1)s0s',
2784        '(a00)','(a00)00s','(a1/20)','(a1/20)00s',],
2785#36
2786    'C m c 21':['(00g)','(00g)s0s','(10g)','(10g)s0s','(a00)','(a00)0s0','(0b0)','(0b0)s00',],
2787    'A 21 m a':['(a00)','(a00)ss0','(a10)','(a10)ss0','(0b0)','(0b0)00s','(00g)','(00g)0s0',],
2788    'B m 21 b':['(0b0)','(0b0)ss0','(1b0)','(1b0)ss0','(a00)','(a00)00s','(00g)','(00g)s00',],
2789    'B b 21 m':['(0b0)','(0b0)0ss','(0b1)','(0b1)ss0','(a00)','(a00)00s','(00g)','(00g)s00',],
2790    'C c m 21':['(00g)','(00g)0ss','(01g)','(01g)0ss','(a00)','(a00)0s0','(0b0)','(0b0)s00',],
2791    'A 21 a m':['(a00)','(a00)s0s','(a01)','(a01)s0s','(0b0)','(0b0)00s','(00g)','(00g)0s0',],
2792#37
2793    'C c c 2':['(00g)','(00g)s0s','(00g)0ss','(10g)','(10g)s0s','(10g)0ss','(01g)','(01g)s0s','(01g)0ss',
2794        '(a00)','(a00)0s0','(0b0)','(0b0)s00',],
2795    'A 2 a a':['(a00)','(a00)ss0','(a00)s0s','(a10)','(a10)ss0','(a10)ss0','(a01)','(a01)ss0','(a01)ss0',
2796        '(0b0)','(0b0)00s','(00g)','(00g)0s0',],
2797    'B b 2 b':['(0b0)','(0b0)0ss','(0b0)ss0','(0b1)','(0b1)0ss','(0b1)ss0','(1b0)','(1b0)0ss','(1b0)ss0',
2798        '(a00)','(a00)00s','(00g)','(00g)s00',],
2799#38
2800    'A m m 2':['(a00)','(a00)0s0','(a10)','(a10)0s0','(00g)','(00g)0s0',
2801        '(00g)ss0','(00g)0ss','(1/20g)','(1/20g)0ss','(0b0)','(0b0)s00','(1/2b0)',],
2802    'B 2 m m':['(0b0)','(0b0)00s','(0b1)','(0b1)00s','(a00)','(a00)00s',
2803        '(a00)0ss','(a00)s0s','(a1/20)','(a1/20)s0s','(00g)','(00g)0s0','(01/2g)',],
2804    'C m 2 m':['(00g)','(00g)s00','(10g)','(10g)s00','(0b0)','(0b0)s00',
2805        '(0b0)s0s','(0b0)ss0','(0b1/2)','(0b1/2)ss0','(a00)','(a00)00s','(a01/2)',],
2806    'A m 2 m':['(a00)','(a00)00s','(a01)','(a01)00s','(0b0)','(0b0)00s',
2807        '(0b0)s0s','(0b0)0ss','(1/2b0)','(1/2b0)0ss','(00g)','(00g)s00','(1/20g)',],
2808    'B m m 2':['(0b0)','(0b0)s00','(0b1)','(0b1)s00','(a00)','(a00)0s0',
2809        '(a00)0ss','(a00)ss0','(01/2g)','(01/2g)s0s','(a00)','(a00)0s0','(a1/20)',],
2810    'C 2 m m':['(00g)','(00g)0s0','(10g)','(10g)0s0','(00g)','(00g)s00',
2811        '(0b0)s0s','(0b0)0ss','(a01/2)','(a01/2)ss0','(0b0)','(0b0)00s','(0b1/2)',],
2812#39
2813    'A b m 2':['(a00)','(a00)0s0','(a01)','(a01)0s0','(00g)','(00g)s0s',
2814        '(00g)ss0','(00g)0ss','(1/20g)','(1/20g)0ss','(0b0)','(0b0)s00','(1/2b0)',],
2815    'B 2 c m':['(0b0)','(0b0)00s','(1b0)','(1b0)00s','(a00)','(a00)ss0',
2816        '(a00)0ss','(a00)s0s','(a1/20)','(a1/20)s0s','(00g)','(00g)0s0','(01/2g)',],
2817    'C m 2 a':['(00g)','(00g)s00','(01g)','(01g)s00','(0b0)','(0b0)0ss',
2818        '(0b0)s0s','(0b0)ss0','(0b1/2)','(0b1/2)ss0','(a00)','(a00)00s','(a01/2)',],
2819    'A c 2 m':['(a00)','(a00)00s','(a10)','(a10)00s','(0b0)','(0b0)ss0',
2820        '(0b0)s0s','(0b0)0ss','(1/2b0)','(1/2b0)0ss','(00g)','(00g)s00','(1/20g)',],
2821    'B m a 2':['(0b0)','(0b0)s00','(0b1)','(0b1)s00','(00g)','(00g)s0s',
2822        '(00g)0ss','(00g)ss0','(01/2g)','(01/2g)ss0','(a00)','(a00)00s','(a1/20)',],
2823    'C 2 m b':['(00g)','(00g)0s0','(10g)','(10g)0s0','(a00)','(a00)0ss',
2824        '(a00)ss0','(a00)s0s','(a01/2)','(a01/2)s0s','(0b0)','(0b0)0s0','(0b1/2)',],
2825#40       
2826    'A m a 2':['(a00)','(a01)','(00g)','(00g)s0s','(00g)ss0','(00g)0ss','(0b0)','(0b0)s00',],
2827    'B 2 m b':['(0b0)','(1b0)','(a00)','(a00)ss0','(a00)0ss','(a00)s0s','(00g)','(00g)0s0',],
2828    'C c 2 m':['(00g)','(01g)','(0b0)','(0b0)0ss','(0b0)s0s','(0b0)ss0','(a00)','(a00)00s',],
2829    'A m 2 a':['(a00)','(a10)','(0b0)','(0b0)ss0','(0b0)s0s','(0b0)0ss','(00g)','(00g)s00',],
2830    'B b m 2':['(0b0)','(0b1)','(00g)','(00g)0ss','(00g)ss0','(00g)s0s','(a00)','(a00)0s0',],
2831    'C 2 c m':['(00g)','(10g)','(a00)','(a00)s0s','(a00)0ss','(a00)ss0','(0b0)','(0b0)00s',],
2832#41
2833    'A b a 2':['(a00)','(a01)','(00g)','(00g)s0s','(00g)ss0','(00g)0ss','(0b0)','(0b0)s00',],
2834    'B 2 c b':['(0b0)','(1b0)','(a00)','(a00)ss0','(a00)0ss','(a00)s0s','(00g)','(00g)0s0',],
2835    'C c 2 a':['(00g)','(01g)','(0b0)','(0b0)0ss','(0b0)s0s','(0b0)ss0','(a00)','(a00)00s',],
2836    'A c 2 a':['(a00)','(a10)','(0b0)','(0b0)ss0','(0b0)s0s','(0b0)0ss','(00g)','(00g)s00',],
2837    'B b a 2':['(0b0)','(0b1)','(00g)','(00g)0ss','(00g)ss0','(00g)s0s','(a00)','(a00)0s0',],
2838    'C 2 c b':['(00g)','(10g)','(a00)','(a00)s0s','(a00)0ss','(a00)ss0','(0b0)','(0b0)00s',],
2839       
2840#42       
2841    'F m m 2':['(00g)','(00g)s0s','(00g)0ss','(00g)ss0','(10g)','(10g)ss0','(10g)s0s',
2842        '(01g)','(01g)ss0','(01g)0ss','(a00)','(a00)0s0','(a01)','(a01)0s0',
2843        '(0b0)','(0b0)s00','(0b1)','(0b1)s00',],       
2844    'F 2 m m':['(a00)','(a00)ss0','(a00)s0s','(a00)0ss','(a10)','(a10)0ss','(a10)ss0',
2845        '(a01)','(a01)0ss','(a01)s0s','(0b0)','(0b0)00s','(1b0)','(1b0)00s',
2846        '(00g)','(00g)0s0','(10g)','(10g)0s0',],
2847    'F m 2 m':['(0b0)','(0b0)0ss','(0b0)ss0','(0b0)s0s','(0b1)','(0b1)s0s','(0b1)0ss',
2848        '(1b0)','(1b0)s0s','(1b0)ss0','(00g)','(00g)s00','(01g)','(01g)s00',
2849        '(a00)','(a00)00s','(a10)','(a10)00s',],       
2850#43       
2851    'F d d 2':['(00g)','(00g)0ss','(00g)s0s','(a00)','(0b0)',],
2852    'F 2 d d':['(a00)','(a00)s0s','(a00)ss0','(00g)','(0b0)',],       
2853    'F d 2 d':['(0b0)','(0b0)0ss','(0b0)ss0','(a00)','(00g)',],
2854#44
2855    'I m m 2':['(00g)','(00g)ss0','(00g)s0s','(00g)0ss','(a00)','(a00)0s0','(0b0)','(0b0)s00',],
2856    'I 2 m m':['(a00)','(00g)0ss','(00g)ss0','(00g)s0s','(0b0)','(0b0)00s','(00g)','(00g)0s0',],
2857    'I m 2 m':['(0b0)','(0b0)s0s','(0b0)0ss','(0b0)ss0','(00g)','(00g)s00','(a00)','(a00)00s',],
2858#45       
2859    'I b a 2':['(00g)','(00g)ss0','(00g)s0s','(00g)0ss','(a00)','(a00)0s0','(0b0)','(0b0)s00',],
2860    'I 2 c b':['(0b0)','(0b0)s0s','(0b0)0ss','(0b0)ss0','(00g)','(00g)s00','(a00)','(a00)00s',],
2861    'I c 2 a':['(a00)','(00g)0ss','(00g)ss0','(00g)s0s','(0b0)','(0b0)00s','(00g)','(00g)0s0',],
2862#46       
2863    'I m a 2':['(a00)','(00g)0ss','(00g)ss0','(00g)s0s','(0b0)','(0b0)00s','(00g)','(00g)0s0',],
2864    'I 2 m b':['(0b0)','(0b0)s0s','(0b0)0ss','(0b0)ss0','(00g)','(00g)s00','(a00)','(a00)00s',],       
2865    'I c 2 m':['(a00)','(00g)0ss','(00g)ss0','(00g)s0s','(0b0)','(0b0)00s','(00g)','(00g)0s0',],
2866    'I m 2 a':['(a00)','(00g)0ss','(00g)ss0','(00g)s0s','(0b0)','(0b0)00s','(00g)','(00g)0s0',],
2867    'I b m 2':['(a00)','(00g)0ss','(00g)ss0','(00g)s0s','(0b0)','(0b0)00s','(00g)','(00g)0s0',],
2868    'I 2 c m':['(0b0)','(0b0)s0s','(0b0)0ss','(0b0)ss0','(00g)','(00g)s00','(a00)','(a00)00s',],
2869#47       
2870    'P m m m':['(00g)','(00g)s00','(00g)0s0','(00g)ss0','(01/2g)','(01/2g)s00','(1/20g)','(1/20g)s00','(1/21/2g)',
2871        '(a00)','(a00)0s0','(a00)00s','(a00)0ss','(a01/2)','(a01/2)0s0','(a1/20)','(a1/20)00s','(a1/21/2)',
2872        '(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s','(1/2b0)','(1/2b0)00s','(0b1/2)','(0b1/2)s00','(1/2b1/2)',],
2873#48 o@i qq0,0qq,q0q ->000
2874    'P n n n':['(00g)','(00g)s00','(00g)0s0','(1/21/2g)',
2875        '(a00)','(a00)0s0','(a00)00s','(a1/21/2)',
2876        '(0b0)','(0b0)s00','(0b0)00s','(1/2b1/2)',],
2877#49       
2878    'P c c m':['(00g)','(00g)s00','(00g)0s0','(01/2g)','(1/20g)','(1/21/2g)',
2879        '(a00)','(a00)0s0','(a00)00s','(a00)0ss','(a1/20)','(a1/20)00s',
2880        '(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s','(1/2b0)','(1/2b0)00s',],       
2881    'P m a a':['(a00)','(a00)0s0','(a00)00s','(a01/2)','(a1/20)','(a1/21/2)',
2882        '(0b0)','(0b0)00s','(0b0)s00','(0b0)s0s','(0b1/2)','(0b1/2)s00',
2883        '(00g)','(00g)0s0','(00g)s00','(00g)ss0','(01/2g)','(01/2g)s00',],       
2884    'P b m b':['(0b0)','(0b0)00s','(0b0)s00','(0b1/2)','(1/2b0)','(1/2b1/2)',
2885        '(00g)','(00g)s00','(00g)0s0','(00g)ss0','(1/20g)','(1/20g)0s0',
2886        '(a00)','(a00)00s','(a00)0s0','(a00)0ss','(a01/2)','(a01/2)0s0',],
2887#50 o@i qq0,0qq,q0q ->000
2888    'P b a n':['(00g)','(00g)s00','(00g)0s0','(00g)ss0','(1/21/2g)',
2889        '(a00)','(a00)0s0','(a01/2)','(0b0)','(0b0)s00','(0b1/2)',],
2890    'P n c b':['(a00)','(a00)0s0','(a00)00s','(a00)0ss','(a1/21/2)',
2891        '(0b0)','(0b0)00s','(1/2b0)','(00g)','(00g)0s0','(1/20g)',],
2892    'P c n a':['(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s','(1/2b1/2)',
2893        '(00g)','(00g)s00','(01/2g)','(a00)','(a00)00s','(a1/20)',],
2894#51       
2895    'P m m a':['(00g)','(00g)s00','(00g)ss0','(00g)0s0','(0b0)','(0b0)s00',
2896        '(0b0)s0s','(0b0)00s','(a00)','(a00)0s0','(01/2g)','(01/2g)s00',
2897        '(0b1/2)','(0b1/2)s00','(a01/2)','(a01/2)0s0','(1/2b0)','(1/2b1/2)',],
2898    'P b m m':['(a00)','(a00)0s0','(a00)0ss','(a00)00s','(00g)','(00g)0s0',
2899        '(00g)ss0','(00g)s00','(0b0)','(0b0)00s','(a01/2)','(a01/2)0s0',
2900        '(1/20g)','(1/20g)0s0','(1/2b0)','(1/2b0)00s','(01/2g)','(1/21/2g)',],
2901    'P m c m':['(0b0)','(0b0)00s','(0b0)s0s','(0b0)s00','(a00)','(a00)00s',
2902        '(a00)0ss','(a00)0s0','(00g)','(00g)s00','(1/2b0)','(1/2b0)00s',
2903        '(a1/20)','(a1/20)00s','(01/2g)','(01/2g)s00','(a01/2)','(a1/21/2)',],
2904    'P m a m':['(0b0)','(0b0)s00','(0b0)s0s','(0b0)00s','(00g)','(00g)s00',
2905        '(00g)ss0','(00g)0s0','(a00)','(a00)00s','(0b1/2)','(0b1/2)s00',
2906        '(01/2g)','(01/2g)s00','(a1/20)','(a1/20)00s','(1/20g)','(1/21/2g)',],
2907    'P m m b':['(00g)','(00g)0s0','(00g)ss0','(00g)s00','(a00)','(a00)0s0',
2908        '(a00)0ss','(a00)00s','(0b0)','(0b0)s00','(a00)','(a00)0s0',
2909        '(a01/2)','(a01/2)0s0','(0b1/2)','(0b1/2)s00','(a1/20)','(a1/21/2)',],
2910    'P c m m':['(a00)','(a00)00s','(a00)0ss','(a00)0s0','(0b0)','(0b0)00s',
2911        '(0b0)s0s','(0b0)s00','(00g)','(00g)0s0','(0b0)','(0b0)00s',
2912        '(1/2b0)','(1/2b0)00s','(1/20g)','(1/20g)0s0','(0b1/2)','(1/2b1/2)',],
2913#52   o@i qq0,0qq,q0q ->000     
2914    'P n n a':['(00g)','(00g)s00','(00g)0s0','(a00)','(a00)00s',
2915        '(0b0)','(0b0)00s','(a1/21/2)','(1/2b1/2)',],
2916    'P b n n':['(a00)','(a00)0s0','(a00)00s','(0b0)','(0b0)s00',
2917        '(00g)','(00g)s00','(1/2b1/2)','(1/21/2g)',],
2918    'P n c n':['(0b0)','(0b0)s00','(0b0)00s','(00g)','(00g)0s0',
2919        '(a00)','(a00)0s0','(1/21/2g)','(a1/21/2)',],
2920    'P n a n':['(0b0)','(0b0)s00','(0b0)00s','(00g)','(00g)0s0',
2921        '(a00)','(a00)0s0','(1/21/2g)','(a1/21/2)',],
2922    'P n n b':['(00g)','(00g)s00','(00g)0s0','(a00)','(a00)00s',
2923        '(0b0)','(0b0)00s','(a1/21/2)','(1/2b1/2)',],
2924    'P c n n':['(a00)','(a00)0s0','(a00)00s','(0b0)','(0b0)s00',
2925        '(00g)','(00g)s00','(1/2b1/2)','(1/21/2g)',],
2926#53       
2927    'P m n a':['(00g)','(00g)s00','(a00)','(a00)00s','(0b0)','(0b0)00s',
2928        '(0b0)s0s','(0b0)s00','(01/2g)','(01/2g)s00','(a1/20)',],
2929    'P b m n':['(a00)','(a00)0s0','(0b0)','(0b0)s00','(00g)','(00g)s00',
2930        '(00g)ss0','(00g)0s0','(a01/2)','(a01/2)0s0','(0b1/2)',],
2931    'P n c m':['(0b0)','(0b0)00s','(00g)','(00g)0s0','(a00)','(a00)0s0',
2932        '(a00)0ss','(a00)00s','(1/2b0)','(1/2b0)00s','(1/20g)',],
2933    'P m a n':['(0b0)','(0b0)s00','(a00)','(a00)0s0','(00g)','(00g)0s0',
2934        '(00g)ss0','(00g)s00','(0b1/2)','(0b1/2)s00','(a01/2)',],
2935    'P n m b':['(00g)','(00g)0s0','(0b0)','(0b0)00s','(a00)','(a00)00s',
2936        '(a00)0ss','(a00)0s0','(1/20g)','(1/20g)0s0','(1/2b0)',],
2937    'P c n m':['(a00)','(a00)00s','(00g)','(00g)s00','(0b0)','(0b0)s00',
2938        '(0b0)s0s','(0b0)00s','(a1/20)','(a1/20)00s','(01/2g)',],
2939#54       
2940    'P c c a':['(00g)','(00g)s00','(0b0)','(0b0)s00','(a00)','(a00)0s0',
2941        '(a00)0ss','(a00)00s','(01/2g)','(1/2b0)',],
2942    'P b a a':['(a00)','(a00)0s0','(00g)','(00g)0s0','(0b0)','(0b0)00s',
2943        '(0b0)s0s','(0b0)s00','(a01/2)','(01/2g)',],
2944    'P b c b':['(0b0)','(0b0)00s','(a00)','(a00)00s','(00g)','(00g)s00',
2945        '(00g)ss0','(00g)0s0','(1/2b0)','(a01/2)',],
2946    'P b a b':['(0b0)','(0b0)s00','(00g)','(00g)s00','(a00)','(a00)00s',
2947        '(a00)0ss','(a00)0s0','(0b1/2)','(1/20g)',],
2948    'P c c b':['(00g)','(00g)0s0','(a00)','(a00)0s0','(0b0)','(0b0)s00',
2949        '(0b0)s0s','(0b0)00s','(1/20g)','(a1/20)',],
2950    'P c a a':['(a00)','(a00)00s','(0b0)','(0b0)00s','(00g)','(00g)0s0',
2951        '(00g)ss0','(00g)s00','(a1/20)','(0b1/2)',],
2952#55       
2953    'P b a m':['(00g)','(00g)s00','(00g)0s0','(00g)ss0',
2954        '(a00)','(a00)00s','(a01/2)','(0b0)','(0b0)00s','(0b1/2)'],
2955    'P m c b':['(a00)','(a00)0s0','(a00)00s','(a00)0ss',
2956        '(0b0)','(0b0)s00','(1/2b0)','(00g)','(00g)s00','(1/20g)'],
2957    'P c m a':['(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s',
2958        '(a00)','(a00)0s0','(a1/20)','(00g)','(00g)0s0','(01/2g)'],
2959#56       
2960    'P c c n':['(00g)','(00g)s00','(00g)0s0','(a00)','(a00)0s0',
2961        '(0b0)','(0b0)s00'],
2962    'P n a a':['(a00)','(a00)0s0','(a00)00s','(0b0)','(0b0)00s',
2963        '(00g)','(00g)0s0'],
2964    'P b n b':['(0b0)','(0b0)s00','(0b0)00s','(a00)','(a00)00s',
2965        '(00g)','(00g)s00'],
2966#57       
2967    'P c a m':['(00g)','(00g)0s0','(a00)','(a00)00s','(0b0)','(0b0)s00',
2968        '(0b0)ss0','(0b0)00s','(01/2g)','(a1/20)','(a1/20)00s',],
2969    'P m a b':['(a00)','(a00)00s','(0b0)','(0b0)s00','(00g)','(00g)0s0',
2970        '(00g)s0s','(00g)s00','(a01/2)','(0b1/2)','(0b1/2)s00',],
2971    'P c m b':['(0b0)','(0b0)s00','(00g)','(00g)0s0','(a00)','(a00)00s',
2972        '(a00)0ss','(a00)0s0','(1/2b0)','(1/20g)','(1/20g)0s0',],
2973    'P b m a':['(0b0)','(0b0)00s','(a00)','(a00)0s0','(00g)','(00g)s00',
2974        '(00g)ss0','(00g)0s0','(0b1/2)','(a01/2)','(a01/2)0s0',],
2975    'P m c a':['(a00)','(a00)0s0','(00g)','(00g)s00','(0b0)','(0b0)00s',
2976        '(0b0)s0s','(0b0)s00','(a1/20)','(01/2g)','(01/2g)s00'],
2977    'P b c m':['(00g)','(00g)s00','(0b0)','(0b0)00s','(a00)','(a00)0s0',
2978        '(a00)0ss','(a00)00s','(1/20g)','(1/2b0)','(1/2b0)00s',],
2979#58       
2980    'P n n m':['(00g)','(00g)s00','(00g)0s0','(a00)',
2981        '(a00)00s','(0b0)','(0b0)00s'],
2982    'P m n n':['(00g)','(00g)s00','(a00)','(a00)0s0',
2983        '(a00)00s','(0b0)','(0b0)s00'],
2984    'P n m n':['(00g)','(00g)0s0','(a00)','(a00)0s0',
2985        '(0b0)','(0b0)s00','(0b0)00s',],
2986#59 o@i
2987    'P m m n':['(00g)','(00g)s00','(00g)0s0','(00g)ss0','(a00)','(a00)0s0',
2988        '(a01/2)','(a01/2)0s0','(0b0)','(0b0)s00','(0b1/2)','(0b1/2)s00',],
2989    'P n m m':['(a00)','(a00)0s0','(a00)00s','(a00)0ss','(00g)','(00g)0s0',
2990        '(1/20g)','(1/20g)0s0','(0b0)','(0b0)00s','(1/2b0)','(1/2b0)00s'],
2991    'P m n m':['(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s','(00g)','(00g)s00',
2992        '(01/2g)','(01/2g)s00','(a00)','(a00)00s','(a1/20)','(a1/20)00s'],
2993#60       
2994    'P b c n':['(00g)','(00g)s00','(00g)0s0','(a00)','(a00)0s0',
2995        '(a00)00s','(0b0)','(0b0)s00','(0b0)00s'],
2996    'P n c a':['(00g)','(00g)s00','(00g)0s0','(a00)','(a00)0s0',
2997        '(a00)00s','(0b0)','(0b0)s00','(0b0)00s'],
2998    'P b n a':['(00g)','(00g)s00','(00g)0s0','(a00)','(a00)0s0',
2999        '(a00)00s','(0b0)','(0b0)s00','(0b0)00s'],
3000    'P c n b':['(00g)','(00g)s00','(00g)0s0','(a00)','(a00)0s0',
3001        '(a00)00s','(0b0)','(0b0)s00','(0b0)00s'],
3002    'P c a n':['(00g)','(00g)s00','(00g)0s0','(a00)','(a00)0s0',
3003        '(a00)00s','(0b0)','(0b0)s00','(0b0)00s'],
3004    'P n a b':['(00g)','(00g)s00','(00g)0s0','(a00)','(a00)0s0',
3005        '(a00)00s','(0b0)','(0b0)s00','(0b0)00s'],
3006#61       
3007    'P b c a':['(00g)','(00g)s00','(00g)0s0','(a00)','(a00)0s0','(a00)00s',
3008        '(0b0)','(0b0)s00','(0b0)00s'],
3009    'P c a b':['(00g)','(00g)s00','(00g)0s0','(a00)','(a00)0s0','(a00)00s',
3010        '(0b0)','(0b0)s00','(0b0)00s'],
3011#62       
3012    'P n m a':['(00g)','(00g)0s0','(a00)','(a00)0s0','(0b0)','(0b0)00s'],
3013    'P b n m':['(00g)','(00g)s00','(a00)','(a00)00s','(0b0)','(0b0)00s'],
3014    'P m c n':['(00g)','(00g)s00','(a00)','(a00)0s0','(0b0)','(0b0)s00'],
3015    'P n a m':['(00g)','(00g)0s0','(a00)','(a00)00s','(0b0)','(0b0)00s'],
3016    'P m n b':['(00g)','(00g)s00','(a00)','(a00)00s','(0b0)','(0b0)s00'],
3017    'P c m n':['(00g)','(00g)0s0','(a00)','(a00)0s0','(0b0)','(0b0)s00'],
3018#63
3019    'C m c m':['(00g)','(00g)s00','(10g)','(10g)s00','(a00)','(a00)00s','(a00)0ss','(a00)0s0','(0b0)','(0b0)00s','(0b0)s0s','(0b0)s00',],
3020    'A m m a':['(a00)','(a00)0s0','(a10)','(a10)0s0','(0b0)','(0b0)s00','(0b0)s0s','(00g)00s','(00g)','(00g)s00','(00g)ss0','(00g)0s0',],
3021    'B b m m':['(0b0)','(0b0)00s','(0b1)','(0b1)00s','(00g)','(00g)0s0','(00g)ss0','(00g)s00','(a00)','(a00)0s0','(a00)0ss','(a00)00s',],
3022    'B m m b':['(0b0)','(0b0)s00','(1b0)','(1b0)s00','(a00)','(a00)0s0','(a00)0ss','(a00)00s','(00g)','(00g)0s0','(00g)ss0','(00g)s00',],
3023    'C c m m':['(00g)','(00g)0s0','(01g)','(01g)0s0','(0b0)','(0b0)00s','(0b0)s0s','(0b0)s00','(a00)','(a00)00s','(a00)0ss','(a00)0s0',],
3024    'A m a m':['(a00)','(a00)00s','(a01)','(a01)00s','(00g)','(00g)s00','(00g)ss0','(00g)0s0','(0b0)','(0b0)s00','(0b0)s0s','(0b0)00s',],
3025#64       
3026    'C m c a':['(00g)','(00g)s00','(10g)','(10g)s00','(0b0)','(0b0)00s','(0b0)s0s','(0b0)s00','(a00)','(a00)00s','(a00)0ss','(a00)0s0',],
3027    'A b m a':['(a00)','(a00)0s0','(a10)','(a10)0s0','(00g)','(00g)s00','(00g)ss0','(00g)0s0','(0b0)','(0b0)s00','(0b0)s0s','(0b0)00s',],
3028    'B b c m':['(0b0)','(0b0)00s','(0b1)','(0b1)00s','(a00)','(a00)0s0','(a00)0ss','(a00)00s','(00g)','(00g)0s0','(00g)ss0','(00g)s00',],
3029    'B m a b':['(0b0)','(0b0)s00','(1b0)','(1b0)s00','(00g)','(00g)0s0','(00g)ss0','(00g)s00','(a00)','(a00)0s0','(a00)0ss','(a00)00s',],
3030    'C c m b':['(00g)','(00g)0s0','(01g)','(01g)0s0','(a00)','(a00)00s','(a00)0ss','(a00)0s0','(0b0)','(0b0)00s','(0b0)s0s','(0b0)s00',],
3031    'A c a m':['(a00)','(a00)00s','(a01)','(a01)00s','(0b0)','(0b0)s00','(0b0)s0s','(0b0)00s','(00g)','(00g)s00','(00g)ss0','(00g)0s0',],
3032#65       
3033    'C m m m':['(00g)','(00g)s00','(00g)ss0','(10g)','(10g)s00','(10g)ss0','(0b0)','(0b0)00s','(0b0)s0s','(0b0)s00','(0b1/2)','(0b1/2)s00',],
3034    'A m m m':['(a00)','(a00)0s0','(a00)0ss','(a10)','(a10)0s0','(a10)0ss','(00g)','(00g)s00','(00g)ss0','(00g)0s0','(1/20g)','(1/20g)0s0',],
3035    'B m m m':['(0b0)','(0b0)00s','(0b0)s0s','(0b1)','(0b1)00s','(0b1)s0s','(a00)','(a00)0s0','(a00)0ss','(a00)00s','(a1/20)','(a1/20)00s',],
3036#66       
3037    'C c c m':['(00g)','(00g)s00','(10g)','(10g)s00','(0b0)','(0b0)00s','(0b0)s0s','(0b0)s00',],
3038    'A m m a':['(a00)','(a00)0s0','(a10)','(a10)0s0','(00g)','(00g)s00','(00g)ss0','(00g)0s0',],
3039    'B b m b':['(0b0)','(0b0)00s','(0b1)','(0b1)00s','(a00)','(a00)0s0','(a00)0ss','(a00)00s',],
3040#67       
3041    'C m m a':['(00g)','(00g)s00','(00g)ss0','(10g)','(10g)s00','(10g)ss0','(a00)','(a00)00s','(a00)0ss','(a00)0s0','(a01/2)','(a01/2)0s0',],
3042    'A b m m':['(a00)','(a00)0s0','(a00)0ss','(a10)','(a10)0s0','(a10)0ss','(0b0)','(0b0)s00','(0b0)s0s','(0b0)00s','(1/2b0)','(1/2b0)00s',],
3043    'B m c m':['(0b0)','(0b0)00s','(0b0)s0s','(0b1)','(0b1)00s','(0b1)s0s','(00g)','(00g)0s0','(00g)ss0','(00g)s00','(01/2g)','(01/2g)s00',],
3044    'B m a m':['(0b0)','(0b0)s00','(0b0)s0s','(1b0)','(1b0)s00','(1b0)s0s','(a00)','(a00)0s0','(a00)0ss','(a00)00s','(a1/20)','(a1/20)00s',],
3045    'C m m b':['(00g)','(00g)0s0','(00g)ss0','(01g)','(01g)0s0','(01g)ss0','(0b0)','(0b0)00s','(0b0)s0s','(0b0)s00','(0b1/2)','(0b1/2)s00',],
3046    'A c m m':['(a00)','(a00)00s','(a00)0ss','(a01)','(a01)00s','(a01)0ss','(00g)','(00g)s00','(00g)ss0','(00g)0s0','(1/20g)','(1/20g)0s0',],
3047#68 o@i
3048    'C c c a':['(00g)','(00g)s00','(10g)','(01g)','(10g)s00','(01g)s00',
3049        '(a00)','(a00)s00','(a00)ss0','(a00)0s0','(0b0)','(0b0)s00','(0b0)ss0','(0b0)0s0'],
3050    'A b a a':['(a00)','(a00)s00','(a10)','(a01)','(a10)s00','(a01)s00',
3051        '(0b0)','(0b0)s00','(0b0)ss0','(0b0)0s0','(00g)','(00g)s00','(00g)ss0','(00g)0s0'],
3052    'B b c b':['(0b0)','(0b0)s00','(0b1)','(1b0)','(0b1)s00','(1b0)s00',
3053        '(00g)','(00g)s00','(00g)ss0','(0b0)0s0','(a00)','(a00)s00','(a00)ss0','(a00)0s0'],
3054    'B b a b':['(0b0)','(0b0)s00','(1b0)','(0b1)','(1b0)s00','(0b1)s00',
3055        '(a00)','(a00)s00','(a00)ss0','(a00)0s0','(00g)','(00g)s00','(00g)ss0','(00g)0s0'],
3056    'C c c b':['(00g)','(00g)ss0','(01g)','(10g)','(01g)s00','(10g)s00',
3057        '(0b0)','(0b0)s00','(0b0)ss0','(0b0)0s0','(a00)','(a00)s00','(a00)ss0','(a00)0s0'],
3058    'A c a a':['(a00)','(a00)ss0','(a01)','(a10)','(a01)s00','(a10)s00',
3059        '(00g)','(00g)s00','(00g)ss0','(00g)0s0','(0b0)','(0b0)s00','(0b0)ss0','(0b0)0s0'],
3060#69       
3061    'F m m m':['(00g)','(00g)s00','(00g)ss0','(a00)','(a00)s00',
3062        '(a00)ss0','(0b0)','(0b0)s00','(0b0)ss0',
3063        '(10g)','(10g)s00','(10g)ss0','(a10)','(a10)0s0',
3064        '(a10)00s','(a10)0ss','(0b1)','(0b1)s00','(0b1)00s','(0b1)s0s',
3065        '(01g)','(01g)s00','(01g)ss0','(a01)','(a01)0s0',
3066        '(a01)00s','(a01)0ss','(1b0)','(1b0)s00','(1b0)00s','(1b0)s0s'],
3067#70 o@i       
3068    'F d d d':['(00g)','(00g)s00','(a00)','(a00)s00','(0b0)','(0b0)s00'],       
3069#71
3070    'I m m m':['(00g)','(00g)s00','(00g)ss0','(a00)','(a00)0s0',
3071        '(a00)ss0','(0b0)','(0b0)s00','(0b0)ss0'],
3072#72       
3073    'I b a m':['(00g)','(00g)s00','(00g)0s0','(00g)ss0','(a00)','(a00)0s0',
3074        '(a00)00s','(a00)0ss','(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s'],
3075    'I m c b':['(a00)','(a00)0s0','(a00)00s','(a00)0ss','(0b0)','(0b0)00s',
3076        '(0b0)s00','(0b0)s0s','(00g)','(00g)0s0','(00g)s00','(00g)ss0'],
3077    'I c m a':['(0b0)','(0b0)00s','(0b0)s00','(0b0)s0s','(00g)','(00g)s00',
3078        '(00g)0s0','(00g)ss0','(a00)','(a00)00s','(a00)0s0','(a00)0ss'],
3079#73       
3080    'I b c a':['(00g)','(00g)s00','(00g)0s0','(00g)ss0','(a00)','(a00)0s0',
3081        '(a00)00s','(a00)0ss','(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s'],
3082    'I c a b':['(00g)','(00g)s00','(00g)0s0','(00g)ss0','(a00)','(a00)0s0',
3083        '(a00)00s','(a00)0ss','(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s'],
3084#74       
3085    'I m m a':['(00g)','(00g)s00','(00g)0s0','(00g)ss0','(a00)','(a00)0s0',
3086        '(a00)00s','(a00)0ss','(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s'],
3087    'I b m m':['(00g)','(00g)s00','(00g)0s0','(00g)ss0','(a00)','(a00)0s0',
3088        '(a00)00s','(a00)0ss','(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s'],
3089    'I m c m':['(00g)','(00g)s00','(00g)0s0','(00g)ss0','(a00)','(a00)0s0',
3090        '(a00)00s','(a00)0ss','(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s'],
3091    'I m a m':['(00g)','(00g)s00','(00g)0s0','(00g)ss0','(a00)','(a00)0s0',
3092        '(a00)00s','(a00)0ss','(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s'],
3093    'I m m b':['(00g)','(00g)s00','(00g)0s0','(00g)ss0','(a00)','(a00)0s0',
3094        '(a00)00s','(a00)0ss','(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s'],
3095    'I c m m':['(00g)','(00g)s00','(00g)0s0','(00g)ss0','(a00)','(a00)0s0',
3096        '(a00)00s','(a00)0ss','(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s'],
3097#tetragonal - done & checked
3098#75
3099    'P 4':['(00g)','(00g)q','(00g)s','(1/21/2g)','(1/21/2g)q',],
3100#76
3101    'P 41':['(00g)','(1/21/2g)',],
3102#77
3103    'P 42':['(00g)','(00g)q','(1/21/2g)','(1/21/2g)q',],
3104#78
3105    'P 43':['(00g)','(1/21/2g)',],
3106#79
3107    'I 4':['(00g)','(00g)q','(00g)s',],
3108#80
3109    'I 41':['(00g)','(00g)q',],
3110#81
3111    'P -4':['(00g)','(1/21/2g)',],
3112#82
3113    'I -4':['(00g)',],
3114#83
3115    'P 4/m':['(00g)','(00g)s0','(1/21/2g)',],
3116#84
3117    'P 42/m':['(00g)','(1/21/2g)',],
3118#85 o@i q0 -> 00
3119    'P 4/n':['(00g)','(00g)s0','(1/21/2g)',], #q0?
3120#86 o@i q0 -> 00
3121    'P 42/n':['(00g)','(1/21/2g)',],      #q0?
3122#87
3123    'I 4/m':['(00g)','(00g)s0',],
3124#88
3125    'I 41/a':['(00g)',],
3126#89
3127    'P 4 2 2':['(00g)','(00g)q00','(00g)s00','(1/21/2g)','(1/21/2g)q00',],
3128#90
3129    'P 4 21 2':['(00g)','(00g)q00','(00g)s00',],
3130#91
3131    'P 41 2 2':['(00g)','(1/21/2g)',],
3132#92
3133    'P 41 21 2':['(00g)',],
3134#93
3135    'P 42 2 2':['(00g)','(00g)q00','(1/21/2g)','(1/21/2g)q00',],
3136#94
3137    'P 42 21 2':['(00g)','(00g)q00',],
3138#95
3139    'P 43 2 2':['(00g)','(1/21/2g)',],
3140#96
3141    'P 43 21 2':['(00g)',],
3142#97
3143    'I 4 2 2':['(00g)','(00g)q00','(00g)s00',],
3144#98
3145    'I 41 2 2':['(00g)','(00g)q00',],
3146#99
3147    'P 4 m m':['(00g)','(00g)ss0','(00g)0ss','(00g)s0s','(1/21/2g)','(1/21/2g)0ss'],
3148#100
3149    'P 4 b m':['(00g)','(00g)ss0','(00g)0ss','(00g)s0s','(1/21/2g)qq0','(1/21/2g)qqs',],
3150#101
3151    'P 42 c m':['(00g)','(00g)0ss','(1/21/2g)','(1/21/2g)0ss',],
3152#102
3153    'P 42 n m':['(00g)','(00g)0ss','(1/21/2g)qq0','(1/21/2g)qqs',],
3154#103
3155    'P 4 c c':['(00g)','(00g)ss0','(1/21/2g)',],
3156#104
3157    'P 4 n c':['(00g)','(00g)ss0','(1/21/2g)qq0',],
3158#105
3159    'P 42 m c':['(00g)','(00g)ss0','(1/21/2g)',],
3160#106
3161    'P 42 b c':['(00g)','(00g)ss0','(1/21/2g)qq0',],
3162#107
3163    'I 4 m m':['(00g)','(00g)ss0','(00g)0ss','(00g)s0s',],
3164#108
3165    'I 4 c m':['(00g)','(00g)ss0','(00g)0ss','(00g)s0s',],
3166#109
3167    'I 41 m d':['(00g)','(00g)ss0',],
3168#110
3169    'I 41 c d':['(00g)','(00g)ss0',],
3170#111
3171    'P -4 2 m':['(00g)','(00g)0ss','(1/21/2g)','(1/21/2g)0ss',],
3172#112
3173    'P -4 2 c':['(00g)','(1/21/2g)',],
3174#113
3175    'P -4 21 m':['(00g)','(00g)0ss',],
3176#114
3177    'P -4 21 c':['(00g)',],
3178#115    00s -> 0ss
3179    'P -4 m 2':['(00g)','(00g)0s0','(1/21/2g)',],
3180#116
3181    'P -4 c 2':['(00g)','(1/21/2g)',],
3182#117    00s -> 0ss
3183    'P -4 b 2':['(00g)','(00g)0s0','(1/21/2g)0q0',],
3184#118
3185    'P -4 n 2':['(00g)','(1/21/2g)0q0',],
3186#119
3187    'I -4 m 2':['(00g)','(00g)0s0',],
3188#120
3189    'I -4 c 2':['(00g)','(00g)0s0',],
3190#121    00s -> 0ss
3191    'I -4 2 m':['(00g)','(00g)0ss',],
3192#122
3193    'I -4 2 d':['(00g)',],
3194#123
3195    'P 4/m m m':['(00g)','(00g)s0s0','(00g)00ss','(00g)s00s',
3196        '(1/21/2g)','(1/21/2g)s0s0','(1/21/2g)00ss','(1/21/2g)s00s',],
3197#124
3198    'P 4/m c c':['(00g)','(00g)s0s0','(1/21/2g)',],
3199#125    o@i q0q0 -> 0000, q0qs -> 00ss
3200    'P 4/n b m':['(00g)','(00g)s0s0','(00g)00ss','(00g)s00s','(1/21/2g)','(1/21/2g)00ss',],
3201#126    o@i q0q0 -> 0000
3202    'P 4/n n c':['(00g)','(00g)s0s0','(1/21/2g)',],
3203#127
3204    'P 4/m b m':['(00g)','(00g)s0s0','(00g)00ss','(00g)s00s',],
3205#128
3206    'P 4/m n c':['(00g)','(00g)s0s0',],
3207#129
3208    'P 4/n m m':['(00g)','(00g)s0s0','(00g)00ss','(00g)s00s',],
3209#130
3210    'P 4/n c c':['(00g)','(00g)s0s0',],
3211#131
3212    'P 42/m m c':['(00g)','(00g)s0s0','(1/21/2g)',],
3213#132
3214    'P 42/m c m':['(00g)','(00g)00ss','(1/21/2g)','(1/21/2g)00ss',],
3215#133    o@i q0q0 -> 0000
3216    'P 42/n b c':['(00g)','(00g)s0s0','(1/21/2g)',],
3217#134    o@i q0q0 -> 0000, q0qs -> 00ss
3218    'P 42/n n m':['(00g)','(00g)00ss','(1/21/2g)','(1/21/2g)00ss',],
3219#135
3220    'P 42/m b c':['(00g)','(00g)s0s0',],
3221#136
3222    'P 42/m n m':['(00g)','(00g)00ss',],
3223#137
3224    'P 42/n m c':['(00g)','(00g)s0s0',],
3225#138
3226    'P 42/n c m':['(00g)','(00g)00ss',],
3227#139
3228    'I 4/m m m':['(00g)','(00g)s0s0','(00g)00ss','(00g)s00s',],
3229#140
3230    'I 4/m c m':['(00g)','(00g)s0s0','(00g)00ss','(00g)s00s',],
3231#141
3232    'I 41/a m d':['(00g)','(00g)s0s0',],
3233#142
3234    'I 41/a c d':['(00g)','(00g)s0s0',],
3235    #trigonal/rhombahedral - done & checked
3236#143
3237    'P 3':['(00g)','(00g)t','(1/31/3g)',],
3238#144
3239    'P 31':['(00g)','(1/31/3g)',],
3240#145
3241    'P 32':['(00g)','(1/31/3g)',],
3242#146
3243    'R 3':['(00g)','(00g)t',],
3244#147
3245    'P -3':['(00g)','(1/31/3g)',],
3246#148
3247    'R -3':['(00g)',],
3248#149
3249    'P 3 1 2':['(00g)','(00g)t00','(1/31/3g)',],
3250#150
3251    'P 3 2 1':['(00g)','(00g)t00',],
3252#151
3253    'P 31 1 2':['(00g)','(1/31/3g)',],
3254#152
3255    'P 31 2 1':['(00g)',],
3256#153
3257    'P 32 1 2':['(00g)','(1/31/3g)',],
3258#154
3259    'P 32 2 1':['(00g)',],
3260#155
3261    'R 3 2':['(00g)','(00g)t0',],
3262#156
3263    'P 3 m 1':['(00g)','(00g)0s0',],
3264#157
3265    'P 3 1 m':['(00g)','(00g)00s','(1/31/3g)','(1/31/3g)00s',],
3266#158
3267    'P 3 c 1':['(00g)',],
3268#159
3269    'P 3 1 c':['(00g)','(1/31/3g)',],
3270#160
3271    'R 3 m':['(00g)','(00g)0s',],
3272#161
3273    'R 3 c':['(00g)',],
3274#162
3275    'P -3 1 m':['(00g)','(00g)00s','(1/31/3g)','(1/31/3g)00s',],
3276#163
3277    'P -3 1 c':['(00g)','(1/31/3g)',],
3278#164
3279    'P -3 m 1':['(00g)','(00g)0s0',],
3280#165
3281    'P -3 c 1':['(00g)',],
3282#166       
3283    'R -3 m':['(00g)','(00g)0s',],
3284#167
3285    'R -3 c':['(00g)',],
3286    #hexagonal - done & checked
3287#168
3288    'P 6':['(00g)','(00g)h','(00g)t','(00g)s',],
3289#169
3290    'P 61':['(00g)',],
3291#170
3292    'P 65':['(00g)',],
3293#171
3294    'P 62':['(00g)','(00g)h',],
3295#172
3296    'P 64':['(00g)','(00g)h',],
3297#173
3298    'P 63':['(00g)','(00g)h',],
3299#174
3300    'P -6':['(00g)',],
3301#175
3302    'P 6/m':['(00g)','(00g)s0',],
3303#176
3304    'P 63/m':['(00g)',],
3305#177
3306    'P 6 2 2':['(00g)','(00g)h00','(00g)t00','(00g)s00',],
3307#178
3308    'P 61 2 2':['(00g)',],
3309#179
3310    'P 65 2 2':['(00g)',],
3311#180
3312    'P 62 2 2':['(00g)','(00g)h00',],
3313#181
3314    'P 64 2 2':['(00g)','(00g)h00',],
3315#182
3316    'P 63 2 2':['(00g)','(00g)h00',],
3317#183
3318    'P 6 m m':['(00g)','(00g)ss0','(00g)0ss','(00g)s0s',],
3319#184
3320    'P 6 c c':['(00g)','(00g)s0s',],
3321#185
3322    'P 63 c m':['(00g)','(00g)0ss',],
3323#186
3324    'P 63 m c':['(00g)','(00g)0ss',],
3325#187
3326    'P -6 m 2':['(00g)','(00g)0s0',],
3327#188
3328    'P -6 c 2':['(00g)',],
3329#189
3330    'P -6 2 m':['(00g)','(00g)00s',],
3331#190
3332    'P -6 2 c':['(00g)',],
3333#191
3334    'P 6/m m m':['(00g)','(00g)s0s0','(00g)00ss','(00g)s00s',],
3335#192
3336    'P 6/m c c':['(00g)','(00g)s00s',],
3337#193
3338    'P 63/m c m':['(00g)','(00g)00ss',],
3339#194
3340    'P 63/m m c':['(00g)','(00g)00ss'],
3341    }
3342
3343#'A few non-standard space groups for test use'
3344nonstandard_sglist = ('P 21 1 1','P 1 21 1','P 1 1 21','R 3 r','R 3 2 h', 
3345                      'R -3 r', 'R 3 2 r','R 3 m h', 'R 3 m r',
3346                      'R 3 c r','R -3 c r','R -3 m r',),
3347
3348#A list of orthorhombic space groups that were renamed in the 2002 Volume A,
3349# along with the pre-2002 name. The e designates a double glide-plane'''
3350sgequiv_2002_orthorhombic= (('A e m 2', 'A b m 2',),
3351                            ('A e a 2', 'A b a 2',),
3352                            ('C m c e', 'C m c a',),
3353                            ('C m m e', 'C m m a',),
3354                            ('C c c e', 'C c c a'),)
3355#Use the space groups types in this order to list the symbols in the
3356#order they are listed in the International Tables, vol. A'''
3357symtypelist = ('triclinic', 'monoclinic', 'orthorhombic', 'tetragonal', 
3358               'trigonal', 'hexagonal', 'cubic')
3359
3360# self-test materials follow. Requires files in directory testinp
3361selftestlist = []
3362'''Defines a list of self-tests'''
3363selftestquiet = True
3364def _ReportTest():
3365    'Report name and doc string of current routine when ``selftestquiet`` is False'
3366    if not selftestquiet:
3367        import inspect
3368        caller = inspect.stack()[1][3]
3369        doc = eval(caller).__doc__
3370        if doc is not None:
3371            print('testing '+__file__+' with '+caller+' ('+doc+')')
3372        else:
3373            print('testing '+__file__()+" with "+caller)
3374def test0():
3375    '''self-test #0: exercise MoveToUnitCell'''
3376    _ReportTest()
3377    msg = "MoveToUnitCell failed"
3378    assert (MoveToUnitCell([1,2,3]) == [0,0,0]).all, msg
3379    assert (MoveToUnitCell([2,-1,-2]) == [0,0,0]).all, msg
3380    assert abs(MoveToUnitCell(np.array([-.1]))[0]-0.9) < 1e-6, msg
3381    assert abs(MoveToUnitCell(np.array([.1]))[0]-0.1) < 1e-6, msg
3382selftestlist.append(test0)
3383
3384def test1():
3385    '''self-test #1: SpcGroup against previous results'''
3386    #'''self-test #1: SpcGroup and SGPrint against previous results'''
3387    _ReportTest()
3388    testdir = ospath.join(ospath.split(ospath.abspath( __file__ ))[0],'testinp')
3389    if ospath.exists(testdir):
3390        if testdir not in sys.path: sys.path.insert(0,testdir)
3391    import spctestinp
3392    def CompareSpcGroup(spc, referr, refdict, reflist): 
3393        'Compare output from GSASIIspc.SpcGroup with results from a previous run'
3394        # if an error is reported, the dictionary can be ignored
3395        msg0 = "CompareSpcGroup failed on space group %s" % spc
3396        result = SpcGroup(spc)
3397        if result[0] == referr and referr > 0: return True
3398        keys = result[1].keys()
3399        #print result[1]['SpGrp']
3400        #msg = msg0 + " in list lengths"
3401        #assert len(keys) == len(refdict.keys()), msg
3402        for key in refdict.keys():
3403            if key == 'SGOps' or  key == 'SGCen':
3404                msg = msg0 + (" in key %s length" % key)
3405                assert len(refdict[key]) == len(result[1][key]), msg
3406                for i in range(len(refdict[key])):
3407                    msg = msg0 + (" in key %s level 0" % key)
3408                    assert np.allclose(result[1][key][i][0],refdict[key][i][0]), msg
3409                    msg = msg0 + (" in key %s level 1" % key)
3410                    assert np.allclose(result[1][key][i][1],refdict[key][i][1]), msg
3411            else:
3412                msg = msg0 + (" in key %s" % key)
3413                assert result[1][key] == refdict[key], msg
3414        msg = msg0 + (" in key %s reflist" % key)
3415        #for (l1,l2) in zip(reflist, SGPrint(result[1])):
3416        #    assert l2.replace('\t','').replace(' ','') == l1.replace(' ',''), 'SGPrint ' +msg
3417        # for now disable SGPrint testing, output has changed
3418        #assert reflist == SGPrint(result[1]), 'SGPrint ' +msg
3419    for spc in spctestinp.SGdat:
3420        CompareSpcGroup(spc, 0, spctestinp.SGdat[spc], spctestinp.SGlist[spc] )
3421selftestlist.append(test1)
3422
3423def test2():
3424    '''self-test #2: SpcGroup against cctbx (sgtbx) computations'''
3425    _ReportTest()
3426    testdir = ospath.join(ospath.split(ospath.abspath( __file__ ))[0],'testinp')
3427    if ospath.exists(testdir):
3428        if testdir not in sys.path: sys.path.insert(0,testdir)
3429    import sgtbxtestinp
3430    def CompareWcctbx(spcname, cctbx_in, debug=0):
3431        'Compare output from GSASIIspc.SpcGroup with results from cctbx.sgtbx'
3432        cctbx = cctbx_in[:] # make copy so we don't delete from the original
3433        spc = (SpcGroup(spcname))[1]
3434        if debug: print spc['SpGrp']
3435        if debug: print spc['SGCen']
3436        latticetype = spcname.strip().upper()[0]
3437        # lattice type of R implies Hexagonal centering", fix the rhombohedral settings
3438        if latticetype == "R" and len(spc['SGCen']) == 1: latticetype = 'P'
3439        assert latticetype == spc['SGLatt'], "Failed: %s does not match Lattice: %s" % (spcname, spc['SGLatt'])
3440        onebar = [1]
3441        if spc['SGInv']: onebar.append(-1)
3442        for (op,off) in spc['SGOps']:
3443            for inv in onebar:
3444                for cen in spc['SGCen']:
3445                    noff = off + cen
3446                    noff = MoveToUnitCell(noff)[0]
3447                    mult = tuple((op*inv).ravel().tolist())
3448                    if debug: print "\n%s: %s + %s" % (spcname,mult,noff)
3449                    for refop in cctbx:
3450                        if debug: print refop
3451                        # check the transform
3452                        if refop[:9] != mult: continue
3453                        if debug: print "mult match"
3454                        # check the translation
3455                        reftrans = list(refop[-3:])
3456                        reftrans = MoveToUnitCell(reftrans)[0]
3457                        if all(abs(noff - reftrans) < 1.e-5):
3458                            cctbx.remove(refop)
3459                            break
3460                    else:
3461                        assert False, "failed on %s:\n\t %s + %s" % (spcname,mult,noff)
3462    for key in sgtbxtestinp.sgtbx:
3463        CompareWcctbx(key, sgtbxtestinp.sgtbx[key])
3464selftestlist.append(test2)
3465
3466def test3(): 
3467    '''self-test #3: exercise SytSym (includes GetOprPtrName, GenAtom, GetKNsym)
3468     for selected space groups against info in IT Volume A '''
3469    _ReportTest()
3470    def ExerciseSiteSym (spc, crdlist):
3471        'compare site symmetries and multiplicities for a specified space group'
3472        msg = "failed on site sym test for %s" % spc
3473        (E,S) = SpcGroup(spc)
3474        assert not E, msg
3475        for t in crdlist:
3476            symb, m = SytSym(t[0],S)
3477            if symb.strip() != t[2].strip() or m != t[1]:
3478                print spc,t[0],m,symb,t[2]
3479            assert m == t[1]
3480            #assert symb.strip() == t[2].strip()
3481
3482    ExerciseSiteSym('p 1',[
3483            ((0.13,0.22,0.31),1,'1'),
3484            ((0,0,0),1,'1'),
3485            ])
3486    ExerciseSiteSym('p -1',[
3487            ((0.13,0.22,0.31),2,'1'),
3488            ((0,0.5,0),1,'-1'),
3489            ])
3490    ExerciseSiteSym('C 2/c',[
3491            ((0.13,0.22,0.31),8,'1'),
3492            ((0.0,.31,0.25),4,'2(y)'),
3493            ((0.25,.25,0.5),4,'-1'),
3494            ((0,0.5,0),4,'-1'),
3495            ])
3496    ExerciseSiteSym('p 2 2 2',[
3497            ((0.13,0.22,0.31),4,'1'),
3498            ((0,0.5,.31),2,'2(z)'),
3499            ((0.5,.31,0.5),2,'2(y)'),
3500            ((.11,0,0),2,'2(x)'),
3501            ((0,0.5,0),1,'222'),
3502            ])
3503    ExerciseSiteSym('p 4/n',[
3504            ((0.13,0.22,0.31),8,'1'),
3505            ((0.25,0.75,.31),4,'2(z)'),
3506            ((0.5,0.5,0.5),4,'-1'),
3507            ((0,0.5,0),4,'-1'),
3508            ((0.25,0.25,.31),2,'4(001)'),
3509            ((0.25,.75,0.5),2,'-4(001)'),
3510            ((0.25,.75,0.0),2,'-4(001)'),
3511            ])
3512    ExerciseSiteSym('p 31 2 1',[
3513            ((0.13,0.22,0.31),6,'1'),
3514            ((0.13,0.0,0.833333333),3,'2(100)'),
3515            ((0.13,0.13,0.),3,'2(110)'),
3516            ])
3517    ExerciseSiteSym('R 3 c',[
3518            ((0.13,0.22,0.31),18,'1'),
3519            ((0.0,0.0,0.31),6,'3'),
3520            ])
3521    ExerciseSiteSym('R 3 c R',[
3522            ((0.13,0.22,0.31),6,'1'),
3523            ((0.31,0.31,0.31),2,'3(111)'),
3524            ])
3525    ExerciseSiteSym('P 63 m c',[
3526            ((0.13,0.22,0.31),12,'1'),
3527            ((0.11,0.22,0.31),6,'m(100)'),
3528            ((0.333333,0.6666667,0.31),2,'3m(100)'),
3529            ((0,0,0.31),2,'3m(100)'),
3530            ])
3531    ExerciseSiteSym('I a -3',[
3532            ((0.13,0.22,0.31),48,'1'),
3533            ((0.11,0,0.25),24,'2(x)'),
3534            ((0.11,0.11,0.11),16,'3(111)'),
3535            ((0,0,0),8,'-3(111)'),
3536            ])
3537selftestlist.append(test3)
3538
3539if __name__ == '__main__':
3540    # run self-tests
3541    selftestquiet = False
3542    for test in selftestlist:
3543        test()
3544    print "OK"
Note: See TracBrowser for help on using the repository browser.