source: trunk/GSASIIspc.py @ 1953

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

code simplification in FillUnitCell? & MoveToUnitCell?
fix error in GenAtom?
  • Property svn:eol-style set to native
  • Property svn:keywords set to Date Author Revision URL Id
File size: 156.0 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-10 14:58:17 +0000 (Mon, 10 Aug 2015) $
12# $Author: vondreele $
13# $Revision: 1953 $
14# $URL: trunk/GSASIIspc.py $
15# $Id: GSASIIspc.py 1953 2015-08-10 14:58:17Z 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: 1953 $")
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 = (xyz+10.)%1.
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    for ic,cen in enumerate(SGData['SGCen']):
1114        C = np.array(cen)
1115        for invers in range(int(SGData['SGInv']+1)):
1116            for io,[M,T] in enumerate(SGData['SGOps']):
1117                idup = ((io+1)+100*ic)*(1-2*invers)
1118                XT = np.inner(M,X)+T
1119                if len(Uij):
1120                    U = Uij2U(Uij)
1121                    U = np.inner(M,np.inner(U,M).T)
1122                    newUij = U2Uij(U)
1123                if invers:
1124                    XT = -XT
1125                XT += C
1126                cell = np.zeros(3)
1127                cellj = np.zeros(3)
1128                if Move:
1129                    newX,cellj = MoveToUnitCell(XT)
1130                else:
1131                    newX = XT
1132                cell += (celli+cellj)
1133                if All:
1134                    if np.allclose(newX,X,atol=0.0002):
1135                        idup = False
1136                else:
1137                    if True in [np.allclose(newX,oldX,atol=0.0002) for oldX in XYZEquiv]:
1138                        idup = False
1139                if All or idup:
1140                    XYZEquiv.append(newX)
1141                    Idup.append(idup)
1142                    Cell.append(cell)
1143                    if len(Uij):
1144                        UijEquiv.append(newUij)                   
1145    if len(Uij):
1146        return zip(XYZEquiv,UijEquiv,Idup,Cell)
1147    else:
1148        return zip(XYZEquiv,Idup,Cell)
1149
1150def GenHKLf(HKL,SGData):
1151    '''
1152    Uses old GSAS Fortran routine genhkl.for
1153
1154    :param HKL:  [h,k,l] must be integral values for genhkl.for to work
1155    :param SGData: space group data obtained from SpcGroup
1156    :returns: iabsnt,mulp,Uniq,phi
1157
1158     *   iabsnt = True if reflection is forbidden by symmetry
1159     *   mulp = reflection multiplicity including Friedel pairs
1160     *   Uniq = numpy array of equivalent hkl in descending order of h,k,l
1161     *   phi = phase offset for each equivalent h,k,l
1162
1163    '''
1164    hklf = list(HKL)+[0,]       #could be numpy array!
1165    Ops = SGData['SGOps']
1166    OpM = np.array([op[0] for op in Ops])
1167    OpT = np.array([op[1] for op in Ops])
1168    Inv = SGData['SGInv']
1169    Cen = np.array([cen for cen in SGData['SGCen']])
1170   
1171    Nuniq,Uniq,iabsnt,mulp = pyspg.genhklpy(hklf,len(Ops),OpM,OpT,SGData['SGInv'],len(Cen),Cen)
1172    h,k,l,f = Uniq
1173    Uniq=np.array(zip(h[:Nuniq],k[:Nuniq],l[:Nuniq]))
1174    phi = f[:Nuniq]
1175   
1176    return iabsnt,mulp,Uniq,phi
1177   
1178def checkSSLaue(HKL,SGData,SSGData):
1179    #Laue check here - Toss HKL if outside unique Laue part
1180    h,k,l,m = HKL
1181    if SGData['SGLaue'] == '2/m':
1182        if SGData['SGUniq'] == 'a':
1183            if 'a' in SSGData['modSymb'] and h == 0 and m < 0:
1184                return False
1185            elif 'b' in SSGData['modSymb'] and k == 0 and l ==0 and m < 0:
1186                return False
1187            else:
1188                return True
1189        elif SGData['SGUniq'] == 'b':
1190            if 'b' in SSGData['modSymb'] and k == 0 and m < 0:
1191                return False
1192            elif 'a' in SSGData['modSymb'] and h == 0 and l ==0 and m < 0:
1193                return False
1194            else:
1195                return True
1196        elif SGData['SGUniq'] == 'c':
1197            if 'g' in SSGData['modSymb'] and l == 0 and m < 0:
1198                return False
1199            elif 'a' in SSGData['modSymb'] and h == 0 and k ==0 and m < 0:
1200                return False
1201            else:
1202                return True
1203    elif SGData['SGLaue'] == 'mmm':
1204        if 'a' in SSGData['modSymb']:
1205            if h == 0 and m < 0:
1206                return False
1207            else:
1208                return True
1209        elif 'b' in SSGData['modSymb']:
1210            if k == 0 and m < 0:
1211                return False
1212            else:
1213                return True
1214        elif 'g' in SSGData['modSymb']:
1215            if l == 0 and m < 0:
1216                return False
1217            else:
1218                return True
1219    else:   #tetragonal, trigonal, hexagonal (& triclinic?)
1220        if l == 0 and m < 0:
1221            return False
1222        else:
1223            return True
1224       
1225   
1226def checkSSextc(HKL,SSGData):
1227    Ops = SSGData['SSGOps']
1228    OpM = np.array([op[0] for op in Ops])
1229    OpT = np.array([op[1] for op in Ops])
1230    HKLS = np.array([HKL,-HKL])     #Freidel's Law
1231    DHKL = np.reshape(np.inner(HKLS,OpM)-HKL,(-1,4))
1232    PHKL = np.reshape(np.inner(HKLS,OpT),(-1,))
1233    for dhkl,phkl in zip(DHKL,PHKL)[1:]:    #skip identity
1234        if dhkl.any():
1235            continue
1236        else:
1237            if phkl%1.:
1238                return False
1239    return True
1240                                 
1241def GetOprPtrName(key):
1242    'Needs a doc string'
1243    OprPtrName = {
1244        '-6643':[   2,' 1bar ', 1],'6479' :[  10,'  2z  ', 2],'-6479':[   9,'  mz  ', 3],
1245        '6481' :[   7,'  my  ', 4],'-6481':[   6,'  2y  ', 5],'6641' :[   4,'  mx  ', 6],
1246        '-6641':[   3,'  2x  ', 7],'6591' :[  28,' m+-0 ', 8],'-6591':[  27,' 2+-0 ', 9],
1247        '6531' :[  25,' m110 ',10],'-6531':[  24,' 2110 ',11],'6537' :[  61,'  4z  ',12],
1248        '-6537':[  62,' -4z  ',13],'975'  :[  68,' 3+++1',14],'6456' :[ 114,'  3z1 ',15],
1249        '-489' :[  73,' 3+-- ',16],'483'  :[  78,' 3-+- ',17],'-969' :[  83,' 3--+ ',18],
1250        '819'  :[  22,' m+0- ',19],'-819' :[  21,' 2+0- ',20],'2431' :[  16,' m0+- ',21],
1251        '-2431':[  15,' 20+- ',22],'-657' :[  19,' m101 ',23],'657'  :[  18,' 2101 ',24],
1252        '1943' :[  48,' -4x  ',25],'-1943':[  47,'  4x  ',26],'-2429':[  13,' m011 ',27],
1253        '2429' :[  12,' 2011 ',28],'639'  :[  55,' -4y  ',29],'-639' :[  54,'  4y  ',30],
1254        '-6484':[ 146,' 2010 ', 4],'6484' :[ 139,' m010 ', 5],'-6668':[ 145,' 2100 ', 6],
1255        '6668' :[ 138,' m100 ', 7],'-6454':[ 148,' 2120 ',18],'6454' :[ 141,' m120 ',19],
1256        '-6638':[ 149,' 2210 ',20],'6638' :[ 142,' m210 ',21],              #search ends here
1257        '2223' :[  68,' 3+++2',39],
1258        '6538' :[ 106,'  6z1 ',40],'-2169':[  83,' 3--+2',41],'2151' :[  73,' 3+--2',42],
1259        '2205' :[  79,'-3-+-2',43],'-2205':[  78,' 3-+-2',44],'489'  :[  74,'-3+--1',45],
1260        '801'  :[  53,'  4y1 ',46],'1945' :[  47,'  4x3 ',47],'-6585':[  62,' -4z3 ',48],
1261        '6585' :[  61,'  4z3 ',49],'6584' :[ 114,'  3z2 ',50],'6666' :[ 106,'  6z5 ',51],
1262        '6643' :[   1,' Iden ',52],'-801' :[  55,' -4y1 ',53],'-1945':[  48,' -4x3 ',54],
1263        '-6666':[ 105,' -6z5 ',55],'-6538':[ 105,' -6z1 ',56],'-2223':[  69,'-3+++2',57],
1264        '-975' :[  69,'-3+++1',58],'-6456':[ 113,' -3z1 ',59],'-483' :[  79,'-3-+-1',60],
1265        '969'  :[  84,'-3--+1',61],'-6584':[ 113,' -3z2 ',62],'2169' :[  84,'-3--+2',63],
1266        '-2151':[  74,'-3+--2',64],'0':[0,' ????',0]
1267        }
1268    return OprPtrName[key]
1269
1270def GetKNsym(key):
1271    'Needs a doc string'
1272    KNsym = {
1273        '0'         :'    1   ','1'         :'   -1   ','64'        :'    2(x)','32'        :'    m(x)',
1274        '97'        :'  2/m(x)','16'        :'    2(y)','8'         :'    m(y)','25'        :'  2/m(y)',
1275        '2'         :'    2(z)','4'         :'    m(z)','7'         :'  2/m(z)','134217728' :'   2(yz)',
1276        '67108864'  :'   m(yz)','201326593' :' 2/m(yz)','2097152'   :'  2(0+-)','1048576'   :'  m(0+-)',
1277        '3145729'   :'2/m(0+-)','8388608'   :'   2(xz)','4194304'   :'   m(xz)','12582913'  :' 2/m(xz)',
1278        '524288'    :'  2(+0-)','262144'    :'  m(+0-)','796433'    :'2/m(+0-)','1024'      :'   2(xy)',
1279        '512'       :'   m(xy)','1537'      :' 2/m(xy)','256'       :'  2(+-0)','128'       :'  m(+-0)',
1280        '385'       :'2/m(+-0)','76'        :'  mm2(x)','52'        :'  mm2(y)','42'        :'  mm2(z)',
1281        '135266336' :' mm2(yz)','69206048'  :'mm2(0+-)','8650760'   :' mm2(xz)','4718600'   :'mm2(+0-)',
1282        '1156'      :' mm2(xy)','772'       :'mm2(+-0)','82'        :'  222   ','136314944' :'  222(x)',
1283        '8912912'   :'  222(y)','1282'      :'  222(z)','127'       :'  mmm   ','204472417' :'  mmm(x)',
1284        '13369369'  :'  mmm(y)','1927'      :'  mmm(z)','33554496'  :'  4(100)','16777280'  :' -4(100)',
1285        '50331745'  :'4/m(100)','169869394' :'422(100)','84934738'  :'-42m 100','101711948' :'4mm(100)',
1286        '254804095' :'4/mmm100','536870928 ':'  4(010)','268435472' :' -4(010)','805306393' :'4/m (10)',
1287        '545783890' :'422(010)','272891986' :'-42m 010','541327412' :'4mm(010)','818675839' :'4/mmm010',
1288        '2050'      :'  4(001)','4098'      :' -4(001)','6151'      :'4/m(001)','3410'      :'422(001)',
1289        '4818'      :'-42m 001','2730'      :'4mm(001)','8191'      :'4/mmm001','8192'      :'  3(111)',
1290        '8193'      :' -3(111)','2629888'   :' 32(111)','1319040'   :' 3m(111)','3940737'   :'-3m(111)',
1291        '32768'     :'  3(+--)','32769'     :' -3(+--)','10519552'  :' 32(+--)','5276160'   :' 3m(+--)',
1292        '15762945'  :'-3m(+--)','65536'     :'  3(-+-)','65537'     :' -3(-+-)','134808576' :' 32(-+-)',
1293        '67437056'  :' 3m(-+-)','202180097' :'-3m(-+-)','131072'    :'  3(--+)','131073'    :' -3(--+)',
1294        '142737664' :' 32(--+)','71434368'  :' 3m(--+)','214040961' :'-3m(--+)','237650'    :'   23   ',
1295        '237695'    :'   m3   ','715894098' :'   432  ','358068946' :'  -43m  ','1073725439':'   m3m  ',
1296        '68157504'  :' mm2d100','4456464'   :' mm2d010','642'       :' mm2d001','153092172' :'-4m2 100',
1297        '277348404' :'-4m2 010','5418'      :'-4m2 001','1075726335':'  6/mmm ','1074414420':'-6m2 100',
1298        '1075070124':'-6m2 120','1075069650':'   6mm  ','1074414890':'   622  ','1073758215':'   6/m  ',
1299        '1073758212':'   -6   ','1073758210':'    6   ','1073759865':'-3m(100)','1075724673':'-3m(120)',
1300        '1073758800':' 3m(100)','1075069056':' 3m(120)','1073759272':' 32(100)','1074413824':' 32(120)',
1301        '1073758209':'   -3   ','1073758208':'    3   ','1074135143':'mmm(100)','1075314719':'mmm(010)',
1302        '1073743751':'mmm(110)','1074004034':' mm2z100','1074790418':' mm2z010','1073742466':' mm2z110',
1303        '1074004004':'mm2(100)','1074790412':'mm2(010)','1073742980':'mm2(110)','1073872964':'mm2(120)',
1304        '1074266132':'mm2(210)','1073742596':'mm2(+-0)','1073872930':'222(100)','1074266122':'222(010)',
1305        '1073743106':'222(110)','1073741831':'2/m(001)','1073741921':'2/m(100)','1073741849':'2/m(010)',
1306        '1073743361':'2/m(110)','1074135041':'2/m(120)','1075314689':'2/m(210)','1073742209':'2/m(+-0)',
1307        '1073741828':' m(001) ','1073741888':' m(100) ','1073741840':' m(010) ','1073742336':' m(110) ',
1308        '1074003968':' m(120) ','1074790400':' m(210) ','1073741952':' m(+-0) ','1073741826':' 2(001) ',
1309        '1073741856':' 2(100) ','1073741832':' 2(010) ','1073742848':' 2(110) ','1073872896':' 2(120) ',
1310        '1074266112':' 2(210) ','1073742080':' 2(+-0) ','1073741825':'   -1   '
1311        }
1312    return KNsym[key]       
1313
1314def GetNXUPQsym(siteSym):
1315    '''       
1316    The codes XUPQ are for lookup of symmetry constraints for position(X), thermal parm(U) & magnetic moments
1317    (P&Q-not used in GSAS-II)
1318    '''
1319    NXUPQsym = {
1320        '    1   ':(28,29,28,28),'   -1   ':( 1,29,28, 0),'    2(x)':(12,18,12,25),'    m(x)':(25,18,12,25),
1321        '  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),
1322        '    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),
1323        '   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),
1324        '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),
1325        '  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),
1326        '   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),
1327        '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),
1328        ' mm2(yz)':(10,13, 0,-1),'mm2(0+-)':(11,13, 0,-1),' mm2(xz)':( 8,12, 0,-1),'mm2(+0-)':( 9,12, 0,-1),
1329        ' mm2(xy)':( 6,11, 0,-1),'mm2(+-0)':( 7,11, 0,-1),'  222   ':( 1,10, 0,-1),'  222(x)':( 1,13, 0,-1),
1330        '  222(y)':( 1,12, 0,-1),'  222(z)':( 1,11, 0,-1),'  mmm   ':( 1,10, 0,-1),'  mmm(x)':( 1,13, 0,-1),
1331        '  mmm(y)':( 1,12, 0,-1),'  mmm(z)':( 1,11, 0,-1),'  4(100)':(12, 4,12, 0),' -4(100)':( 1, 4,12, 0),
1332        '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),
1333        '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),
1334        '422(010)':( 1, 3, 0,-1),'-42m 010':( 1, 3, 0,-1),'4mm(010)':(13, 3, 0,-1),'4/mmm010':(1, 3, 0,-1,),
1335        '  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),
1336        '-42m 001':( 1, 2, 0,-1),'4mm(001)':(14, 2, 0,-1),'4/mmm001':( 1, 2, 0,-1),'  3(111)':( 2, 5, 2, 0),
1337        ' -3(111)':( 1, 5, 2, 0),' 32(111)':( 1, 5, 0, 2),' 3m(111)':( 2, 5, 0, 2),'-3m(111)':( 1, 5, 0,-1),
1338        '  3(+--)':( 5, 8, 5, 0),' -3(+--)':( 1, 8, 5, 0),' 32(+--)':( 1, 8, 0, 5),' 3m(+--)':( 5, 8, 0, 5),
1339        '-3m(+--)':( 1, 8, 0,-1),'  3(-+-)':( 4, 7, 4, 0),' -3(-+-)':( 1, 7, 4, 0),' 32(-+-)':( 1, 7, 0, 4),
1340        ' 3m(-+-)':( 4, 7, 0, 4),'-3m(-+-)':( 1, 7, 0,-1),'  3(--+)':( 3, 6, 3, 0),' -3(--+)':( 1, 6, 3, 0),
1341        ' 32(--+)':( 1, 6, 0, 3),' 3m(--+)':( 3, 6, 0, 3),'-3m(--+)':( 1, 6, 0,-1),'   23   ':( 1, 1, 0, 0),
1342        '   m3   ':( 1, 1, 0, 0),'   432  ':( 1, 1, 0, 0),'  -43m  ':( 1, 1, 0, 0),'   m3m  ':( 1, 1, 0, 0),
1343        ' mm2d100':(12,13, 0,-1),' mm2d010':(13,12, 0,-1),' mm2d001':(14,11, 0,-1),'-4m2 100':( 1, 4, 0,-1),
1344        '-4m2 010':( 1, 3, 0,-1),'-4m2 001':( 1, 2, 0,-1),'  6/mmm ':( 1, 9, 0,-1),'-6m2 100':( 1, 9, 0,-1),
1345        '-6m2 120':( 1, 9, 0,-1),'   6mm  ':(14, 9, 0,-1),'   622  ':( 1, 9, 0,-1),'   6/m  ':( 1, 9,14,-1),
1346        '   -6   ':( 1, 9,14, 0),'    6   ':(14, 9,14, 0),'-3m(100)':( 1, 9, 0,-1),'-3m(120)':( 1, 9, 0,-1),
1347        ' 3m(100)':(14, 9, 0,14),' 3m(120)':(14, 9, 0,14),' 32(100)':( 1, 9, 0,14),' 32(120)':( 1, 9, 0,14),
1348        '   -3   ':( 1, 9,14, 0),'    3   ':(14, 9,14, 0),'mmm(100)':( 1,14, 0,-1),'mmm(010)':( 1,15, 0,-1),
1349        'mmm(110)':( 1,11, 0,-1),' mm2z100':(14,14, 0,-1),' mm2z010':(14,15, 0,-1),' mm2z110':(14,11, 0,-1),
1350        'mm2(100)':(12,14, 0,-1),'mm2(010)':(13,15, 0,-1),'mm2(110)':( 6,11, 0,-1),'mm2(120)':(15,14, 0,-1),
1351        'mm2(210)':(16,15, 0,-1),'mm2(+-0)':( 7,11, 0,-1),'222(100)':( 1,14, 0,-1),'222(010)':( 1,15, 0,-1),
1352        '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),
1353        '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),
1354        ' m(001) ':(23,16,14,23),' m(100) ':(26,25,12,26),' m(010) ':(27,28,13,27),' m(110) ':(18,19, 6,18),
1355        ' m(120) ':(24,27,15,24),' m(210) ':(25,26,16,25),' m(+-0) ':(17,20, 7,17),' 2(001) ':(14,16,14,23),
1356        ' 2(100) ':(12,25,12,26),' 2(010) ':(13,28,13,27),' 2(110) ':( 6,19, 6,18),' 2(120) ':(15,27,15,24),
1357        ' 2(210) ':(16,26,16,25),' 2(+-0) ':( 7,20, 7,17),'   -1   ':( 1,29,28, 0)
1358        }
1359    return NXUPQsym[siteSym]
1360
1361def GetCSxinel(siteSym): 
1362    'Needs a doc string'
1363    CSxinel = [[],                         # 0th empty - indices are Fortran style
1364        [[0,0,0],[ 0.0, 0.0, 0.0]],      #1  0  0  0
1365        [[1,1,1],[ 1.0, 1.0, 1.0]],      #2  X  X  X
1366        [[1,1,1],[ 1.0, 1.0,-1.0]],      #3  X  X -X
1367        [[1,1,1],[ 1.0,-1.0, 1.0]],      #4  X -X  X
1368        [[1,1,1],[ 1.0,-1.0,-1.0]],      #5 -X  X  X
1369        [[1,1,0],[ 1.0, 1.0, 0.0]],      #6  X  X  0
1370        [[1,1,0],[ 1.0,-1.0, 0.0]],      #7  X -X  0
1371        [[1,0,1],[ 1.0, 0.0, 1.0]],      #8  X  0  X
1372        [[1,0,1],[ 1.0, 0.0,-1.0]],      #9  X  0 -X
1373        [[0,1,1],[ 0.0, 1.0, 1.0]],      #10  0  Y  Y
1374        [[0,1,1],[ 0.0, 1.0,-1.0]],      #11 0  Y -Y
1375        [[1,0,0],[ 1.0, 0.0, 0.0]],      #12  X  0  0
1376        [[0,1,0],[ 0.0, 1.0, 0.0]],      #13  0  Y  0
1377        [[0,0,1],[ 0.0, 0.0, 1.0]],      #14  0  0  Z
1378        [[1,1,0],[ 1.0, 2.0, 0.0]],      #15  X 2X  0
1379        [[1,1,0],[ 2.0, 1.0, 0.0]],      #16 2X  X  0
1380        [[1,1,2],[ 1.0, 1.0, 1.0]],      #17  X  X  Z
1381        [[1,1,2],[ 1.0,-1.0, 1.0]],      #18  X -X  Z
1382        [[1,2,1],[ 1.0, 1.0, 1.0]],      #19  X  Y  X
1383        [[1,2,1],[ 1.0, 1.0,-1.0]],      #20  X  Y -X
1384        [[1,2,2],[ 1.0, 1.0, 1.0]],      #21  X  Y  Y
1385        [[1,2,2],[ 1.0, 1.0,-1.0]],      #22  X  Y -Y
1386        [[1,2,0],[ 1.0, 1.0, 0.0]],      #23  X  Y  0
1387        [[1,0,2],[ 1.0, 0.0, 1.0]],      #24  X  0  Z
1388        [[0,1,2],[ 0.0, 1.0, 1.0]],      #25  0  Y  Z
1389        [[1,1,2],[ 1.0, 2.0, 1.0]],      #26  X 2X  Z
1390        [[1,1,2],[ 2.0, 1.0, 1.0]],      #27 2X  X  Z
1391        [[1,2,3],[ 1.0, 1.0, 1.0]],      #28  X  Y  Z
1392        ]
1393    indx = GetNXUPQsym(siteSym)
1394    return CSxinel[indx[0]]
1395   
1396def GetCSuinel(siteSym):
1397    "returns Uij terms, multipliers, GUI flags & Uiso2Uij multipliers"
1398    CSuinel = [[],                                             # 0th empty - indices are Fortran style
1399        [[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
1400        [[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
1401        [[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
1402        [[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
1403        [[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
1404        [[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
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]],    #7  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]],    #8  A  A  A  D  D -D
1407        [[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
1408        [[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
1409        [[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
1410        [[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
1411        [[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
1412        [[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
1413        [[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
1414        [[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
1415        [[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
1416        [[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
1417        [[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
1418        [[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
1419        [[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
1420        [[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
1421        [[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
1422        [[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
1423        [[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
1424        [[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
1425        [[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
1426        [[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
1427        [[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
1428        ]
1429    indx = GetNXUPQsym(siteSym)
1430    return CSuinel[indx[1]]
1431   
1432#def GetSSfxuinel(waveType,nH,XYZ,SGData,SSGData,debug=False):
1433#   
1434#    def fracCrenel(tau,Toff,Twid):
1435#        Tau = (tau-Toff)%1.
1436#        A = np.where(Tau<Twid,1.,0.)
1437#        return A
1438#       
1439#    def fracFourier(tau,nH,fsin,fcos):
1440#        SA = np.sin(2.*nH*np.pi*tau)
1441#        CB = np.cos(2.*nH*np.pi*tau)
1442#        A = SA[np.newaxis,np.newaxis,:]*fsin[:,:,np.newaxis]
1443#        B = CB[np.newaxis,np.newaxis,:]*fcos[:,:,np.newaxis]
1444#        return A+B
1445#       
1446#    def posFourier(tau,nH,psin,pcos):
1447#        SA = np.sin(2*nH*np.pi*tau)
1448#        CB = np.cos(2*nH*np.pi*tau)
1449#        A = SA[np.newaxis,np.newaxis,:]*psin[:,:,np.newaxis]
1450#        B = CB[np.newaxis,np.newaxis,:]*pcos[:,:,np.newaxis]
1451#        return A+B   
1452#
1453#    def posSawtooth(tau,Toff,slopes):
1454#        Tau = (tau-Toff[:,np.newaxis])%1.
1455#        A = slopes[:,:,np.newaxis]*Tau
1456#        return A
1457#   
1458#    def posZigZag(tau,Toff,slopes):
1459#        Tau = (tau-Toff[:,np.newaxis])%1.
1460#        A = np.where(Tau <= 0.5,slopes[:,:,np.newaxis]*Tau,slopes[:,:,np.newaxis]*(1.-Tau))
1461#        return A
1462#       
1463#    print 'super space group: ',SSGData['SSpGrp']
1464#    CSI = {'Sfrac':[[[1,0],[2,0]],[[1.,0.],[1.,0.]]],
1465#        'Spos':[[[1,0,0],[2,0,0],[3,0,0], [4,0,0],[5,0,0],[6,0,0]],
1466#            [[1.,0.,0.],[1.,0.,0.],[1.,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.]]],    #sin & cos
1467#        'Sadp':[[[1,0,0],[2,0,0],[3,0,0],[4,0,0],[5,0,0],[6,0,0],
1468#            [7,0,0],[8,0,0],[9,0,0],[10,0,0],[11,0,0],[12,0,0]],
1469#            [[1.,0.,0.],[1.,0.,0.],[1.,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.],
1470#            [1.,0.,0.],[1.,0.,0.],[1.,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.]]],
1471#        'Smag':[[[1,0,0],[2,0,0],[3,0,0], [4,0,0],[5,0,0],[6,0,0]],
1472#            [[1.,0.,0.],[1.,0.,0.],[1.,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.]]],}
1473#    xyz = np.array(XYZ)%1.
1474#    xyzt = np.array(XYZ+[0,])%1.
1475#    SGOps = copy.deepcopy(SGData['SGOps'])
1476#    siteSym = SytSym(XYZ,SGData)[0].strip()
1477#    print 'siteSym: ',siteSym
1478#    if siteSym == '1':   #"1" site symmetry
1479#        if debug:
1480#            return CSI,None,None,None,None
1481#        else:
1482#            return CSI
1483#    elif siteSym == '-1':   #"-1" site symmetry
1484#        CSI['Sfrac'][0] = [[1,0],[0,0]]
1485#        CSI['Spos'][0] = [[1,0,0],[2,0,0],[3,0,0], [0,0,0],[0,0,0],[0,0,0]]
1486#        CSI['Sadp'][0] = [[0,0,0],[0,0,0],[0,0,0],[0,0,0],[0,0,0],[0,0,0],
1487#        [1,0,0],[2,0,0],[3,0,0],[4,0,0],[5,0,0],[6,0,0]]
1488#        if debug:
1489#            return CSI,None,None,None,None
1490#        else:
1491#            return CSI
1492#    SSGOps = copy.deepcopy(SSGData['SSGOps'])
1493#    #expand ops to include inversions if any
1494#    if SGData['SGInv']:
1495#        for op,sop in zip(SGData['SGOps'],SSGData['SSGOps']):
1496#            SGOps.append([-op[0],-op[1]%1.])
1497#            SSGOps.append([-sop[0],-sop[1]%1.])
1498#    #build set of sym ops around special poasition       
1499#    SSop = []
1500#    Sop = []
1501#    Sdtau = []
1502#    for iop,Op in enumerate(SGOps):         
1503#        nxyz = (np.inner(Op[0],xyz)+Op[1])%1.
1504#        if np.allclose(xyz,nxyz,1.e-4) and iop and MT2text(Op).replace(' ','') != '-X,-Y,-Z':
1505#            SSop.append(SSGOps[iop])
1506#            Sop.append(SGOps[iop])
1507#            ssopinv = nl.inv(SSGOps[iop][0])
1508#            mst = ssopinv[3][:3]
1509#            epsinv = ssopinv[3][3]
1510#            Sdtau.append(np.sum(mst*(XYZ-SGOps[iop][1])-epsinv*SSGOps[iop][1][3]))
1511#    SdIndx = np.argsort(np.array(Sdtau))     # just to do in sensible order
1512#    OpText =  [MT2text(s).replace(' ','') for s in Sop]         #debug?
1513#    SSOpText = [SSMT2text(ss).replace(' ','') for ss in SSop]   #debug?
1514#    print 'special pos super operators: ',SSOpText
1515#    #setup displacement arrays
1516#    tau = np.linspace(0,1,49,True)
1517#    delt2 = np.eye(2)*0.001
1518#    delt4 = np.eye(4)*0.001
1519#    delt6 = np.eye(6)*0.001
1520#    delt12 = np.eye(12)*0.0001
1521#    #make modulation arrays - one parameter at a time
1522#    #site fractions
1523#    CSI['Sfrac'] = [np.zeros((2),dtype='i'),np.ones(2)]
1524#    if 'Crenel' in waveType:
1525#        dF = fracCrenel(tau,delt2[:1],delt2[1:]).squeeze()
1526#    else:
1527#        dF = fracFourier(tau,nH,delt2[:1],delt2[1:]).squeeze()
1528#    dFT = np.zeros_like(dF)
1529#    #positions       
1530#    if 'Fourier' in waveType:
1531#        dX = posFourier(tau,nH,delt6[:3],delt6[3:]) #+np.array(XYZ)[:,np.newaxis,np.newaxis]
1532#          #3x6x12 modulated position array (X,Spos,tau)& force positive
1533#        CSI['Spos'] = [np.zeros((6,3),dtype='i'),np.zeros((6,3))]
1534#    elif waveType == 'Sawtooth':
1535#        CSI['Spos'] = [np.array([[1,],[2,],[3,],[4,]]),np.array([[1.0,],[1.0,],[1.0,],[1.0,]])]
1536#    elif waveType == 'ZigZag':
1537#        CSI['Spos'] = [np.array([[1,],[2,],[3,],[4,]]),np.array([[1.0,],[1.0,],[1.0,],[1.0,]])]
1538#    #anisotropic thermal motion
1539#    dU = posFourier(tau,nH,delt12[:6],delt12[6:])                  #Uij modulations - 6x12x12 array
1540#    CSI['Sadp'] = [np.zeros((12,3),dtype='i'),np.zeros((12,3))]
1541#       
1542#    FSC = np.ones(2,dtype='i')
1543#    VFSC = np.ones(2)
1544#    XSC = np.ones(6,dtype='i')
1545#    USC = np.ones(12,dtype='i')
1546#    dFTP = []
1547#    dXTP = []
1548#    dUTP = []
1549#    for i in SdIndx:
1550#        sop = Sop[i]
1551#        ssop = SSop[i]
1552#        fsc = np.ones(2,dtype='i')
1553#        xsc = np.ones(6,dtype='i')
1554#        ssopinv = nl.inv(ssop[0])
1555#        mst = ssopinv[3][:3]
1556#        epsinv = ssopinv[3][3]
1557#        sdet = nl.det(sop[0])
1558#        ssdet = nl.det(ssop[0])
1559#        dtau = mst*(XYZ-sop[1])-epsinv*ssop[1][3]
1560#        dT = 1.0
1561#        if np.any(dtau%.5):
1562#            dT = np.tan(np.pi*np.sum(dtau%.5))
1563#        tauT = np.inner(mst,XYZ-sop[1])+epsinv*(tau-ssop[1][3])
1564#        if waveType == 'Fourier':
1565#            dXT = posFourier(np.sort(tauT),nH,delt6[:3],delt6[3:])   #+np.array(XYZ)[:,np.newaxis,np.newaxis]
1566#        elif waveType == 'Sawtooth':
1567#            dXT = posSawtooth(tauT,delt4[0],delt4[1:])+np.array(XYZ)[:,np.newaxis,np.newaxis]
1568#        elif waveType == 'ZigZag':
1569#            dXT = posZigZag(tauT,delt4[0],delt4[1:])+np.array(XYZ)[:,np.newaxis,np.newaxis]           
1570#        dXT = np.inner(sop[0],dXT.T)
1571#        dXT = np.swapaxes(dXT,1,2)
1572#        dXT[:,:3,:] *= ssdet
1573#        dXTP.append(dXT)
1574#        if waveType == 'Fourier':
1575#            if np.any(dtau%.5) and ('1/2' in SSGData['modSymb'] or '1' in SSGData['modSymb']):
1576#                CSI['Spos'] = [[[1,0,0],[2,0,0],[3,0,0], [1,0,0],[2,0,0],[3,0,0]],
1577#                    [[1.,0.,0.],[1.,0.,0.],[1.,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.]]]                   
1578#                if '(x)' in siteSym:
1579#                    CSI['Spos'][1][3:] = [1./dT,0.,0.],[-dT,0.,0.],[-dT,0.,0.]
1580#                    if 'm' in siteSym and len(SdIndx) == 1:
1581#                        CSI['Spos'][1][3:] = [-dT,0.,0.],[1./dT,0.,0.],[1./dT,0.,0.]
1582#                elif '(y)' in siteSym:
1583#                    CSI['Spos'][1][3:] = [-dT,0.,0.],[1./dT,0.,0.],[-dT,0.,0.]
1584#                    if 'm' in siteSym and len(SdIndx) == 1:
1585#                        CSI['Spos'][1][3:] = [1./dT,0.,0.],[-dT,0.,0.],[1./dT,0.,0.]
1586#                elif '(z)' in siteSym:
1587#                    CSI['Spos'][1][3:] = [-dT,0.,0.],[-dT,0.,0.],[1./dT,0.,0.]
1588#                    if 'm' in siteSym and len(SdIndx) == 1:
1589#                        CSI['Spos'][1][3:] = [1./dT,0.,0.],[1./dT,0.,0.],[-dT,0.,0.]
1590#                for i in range(3):
1591#                    if not XSC[i]:
1592#                        CSI['Spos'][0][i] = [0,0,0]
1593#                        CSI['Spos'][1][i] = [0.,0.,0.]
1594#                        CSI['Spos'][0][i+3] = [0,0,0]
1595#                        CSI['Spos'][1][i+3] = [0.,0.,0.]
1596#            else:
1597#                for i in range(3):
1598#                    if np.allclose(dX[i,i,:],dXT[i,i,:]*sdet):
1599#                        xsc[i] = 1
1600#                    else:
1601#                        xsc[i] = 0
1602#                    if np.allclose(dX[i,i+3,:],dXT[i,i+3,:]):
1603#                        xsc[i+3] = 1
1604#                    else:
1605#                        xsc[i+3] = 0
1606#            XSC &= xsc
1607#           
1608#        fsc = np.ones(2,dtype='i')
1609#        if 'Crenel' in waveType:
1610#            dFT = fracCrenel(tauT,delt2[:1],delt2[1:]).squeeze()
1611#            fsc = [1,1]
1612#        else:
1613#            dFT = fracFourier(tauT,nH,delt2[:1],delt2[1:]).squeeze()
1614#            dFT = nl.det(sop[0])*dFT
1615#            dFT = dFT[:,np.argsort(tauT)]
1616#            dFT[0] *= ssdet
1617#            dFT[1] *= sdet
1618#            dFTP.append(dFT)
1619#       
1620#            if np.any(dtau%.5) and ('1/2' in SSGData['modSymb'] or '1' in SSGData['modSymb']):
1621#                fsc = [1,1]
1622#                CSI['Sfrac'] = [[[1,0],[1,0]],[[1.,0.],[1/dT,0.]]]
1623#                for i in range(2):
1624#                    if not FSC[i]:
1625#                        CSI['Sfrac'][0][i] = [0,0]
1626#                        CSI['Sfrac'][1][i] = [0.,0.]
1627#            else:
1628#                for i in range(2):
1629#                    if np.allclose(dF[i,:],dFT[i,:],atol=1.e-6):
1630#                        fsc[i] = 1
1631#                    else:
1632#                        fsc[i] = 0
1633#        FSC &= fsc
1634#           
1635#        usc = np.ones(12,dtype='i')
1636#        dUT = posFourier(tauT,nH,delt12[:6],delt12[6:])                  #Uij modulations - 6x12x49 array
1637#        dUijT = np.rollaxis(np.rollaxis(np.array(Uij2U(dUT)),3),3)    #convert dUT to 12x49x3x3
1638#        dUijT = np.rollaxis(np.inner(np.inner(sop[0],dUijT),sop[0].T),3)
1639#        dUT = np.array(U2Uij(dUijT))
1640#        dUT = dUT[:,:,np.argsort(tauT)]
1641#        dUTP.append(dUT)
1642#        if np.any(dtau%.5) and ('1/2' in SSGData['modSymb'] or '1' in SSGData['modSymb']):
1643#            CSI['Sadp'] = [[[1,0,0],[2,0,0],[3,0,0],[4,0,0],[5,0,0],[6,0,0],
1644#            [1,0,0],[2,0,0],[3,0,0],[4,0,0],[5,0,0],[6,0,0]],
1645#            [[1.,0.,0.],[1.,0.,0.],[1.,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.],
1646#            [1./dT,0.,0.],[1./dT,0.,0.],[1./dT,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.]]]
1647#            if '(x)' in siteSym:
1648#                CSI['Sadp'][1][9:] = [-dT,0.,0.],[-dT,0.,0.],[1./dT,0.,0.]
1649#            elif '(y)' in siteSym:
1650#                CSI['Sadp'][1][9:] = [-dT,0.,0.],[1./dT,0.,0.],[-dT,0.,0.]
1651#            elif '(z)' in siteSym:
1652#                CSI['Sadp'][1][9:] = [1./dT,0.,0.],[-dT,0.,0.],[-dT,0.,0.]
1653#            for i in range(6):
1654#                if not USC[i]:
1655#                    CSI['Sadp'][0][i] = [0,0,0]
1656#                    CSI['Sadp'][1][i] = [0.,0.,0.]
1657#                    CSI['Sadp'][0][i+6] = [0,0,0]
1658#                    CSI['Sadp'][1][i+6] = [0.,0.,0.]
1659#        else:
1660#                       
1661#            for i in range(6):
1662#                if np.allclose(dU[i,i,:],dUT[i,i,:]*sdet):
1663#                    usc[i] = 1
1664#                else:
1665#                    usc[i] = 0
1666#                if np.allclose(dU[i,i+6,:],dUT[i,i+6,:]):
1667#                    usc[i+6] = 1
1668#                else:
1669#                    usc[i+6] = 0
1670#            if '4/m' in siteSym and np.any(dUT[0,1,:]):
1671#                CSI['Sadp'][0][6:8] = [[12,0,0],[12,0,0]]
1672#                if ssop[1][3]:
1673#                    CSI['Sadp'][1][6:8] = [[1.,0.,0.],[-1.,0.,0.]]
1674#                    usc[9] = 1
1675#                else:
1676#                    CSI['Sadp'][1][6:8] = [[1.,0.,0.],[1.,0.,0.]]
1677#                    usc[9] = 0
1678#            elif '4' in siteSym and np.any(dUT[0,1,:]):
1679#                CSI['Sadp'][0][6:8] = [[12,0,0],[12,0,0]]
1680#                CSI['Sadp'][0][:2] = [[11,0,0],[11,0,0]]
1681#                if ssop[1][3]:
1682#                    CSI['Sadp'][1][:2] = [[1.,0.,0.],[-1.,0.,0.]]
1683#                    CSI['Sadp'][1][6:8] = [[1.,0.,0.],[-1.,0.,0.]]
1684#                    usc[3] = 1
1685#                    usc[9] = 1
1686#                else:
1687#                    CSI['Sadp'][1][:2] = [[1.,0.,0.],[1.,0.,0.]]
1688#                    CSI['Sadp'][1][6:8] = [[1.,0.,0.],[1.,0.,0.]]
1689#                    usc[3] = 0               
1690#                    usc[9] = 0
1691#            print SSMT2text(ssop).replace(' ',''),sdet,ssdet,epsinv,usc
1692#        USC &= usc
1693#    print dtau
1694#    print FSC
1695#    print XSC
1696#    print USC
1697#    if not np.any(dtau%.5):
1698#        n = -1
1699#        for i,U in enumerate(USC):
1700#            if U:
1701#                n += 1
1702#                CSI['Sadp'][0][i][0] = n+1
1703#                CSI['Sadp'][1][i][0] = 1.0
1704#        if waveType == 'Fourier':
1705#            n = -1
1706#            for i,X in enumerate(XSC):
1707#                if X:
1708#                    n += 1
1709#                    CSI['Spos'][0][i][0] = n+1
1710#                    CSI['Spos'][1][i][0] = 1.0
1711#        n = -1
1712#        for i,F in enumerate(FSC):
1713#            if F:
1714#                n += 1
1715#                CSI['Sfrac'][0][i] = n+1
1716#                CSI['Sfrac'][1][i] = 1.0
1717#            else:
1718#                CSI['Sfrac'][0][i] = 0
1719#                CSI['Sfrac'][1][i] = 0.           
1720#    if debug:
1721#        return CSI,[tau,tauT],[dF,dFTP],[dX,dXTP],[dU,dUTP]
1722#    else:
1723#        return CSI
1724#
1725def GetSSfxuinel(waveType,nH,XYZ,SGData,SSGData,debug=False):
1726   
1727    def orderParms(CSI):
1728        parms = [0,]
1729        for csi in CSI:
1730            for i in [0,1,2]:
1731                if csi[i] not in parms:
1732                    parms.append(csi[i])
1733        for csi in CSI:
1734            for i in [0,1,2]:
1735                csi[i] = parms.index(csi[i])
1736        return CSI
1737   
1738    def fracCrenel(tau,Toff,Twid):
1739        Tau = (tau-Toff[:,np.newaxis])%1.
1740        A = np.where(Tau<Twid[:,np.newaxis],1.,0.)
1741        return A
1742       
1743    def fracFourier(tau,nH,fsin,fcos):
1744        SA = np.sin(2.*nH*np.pi*tau)
1745        CB = np.cos(2.*nH*np.pi*tau)
1746        A = SA[np.newaxis,np.newaxis,:]*fsin[:,:,np.newaxis]
1747        B = CB[np.newaxis,np.newaxis,:]*fcos[:,:,np.newaxis]
1748        return A+B
1749       
1750    def posFourier(tau,nH,psin,pcos):
1751        SA = np.sin(2*nH*np.pi*tau)
1752        CB = np.cos(2*nH*np.pi*tau)
1753        A = SA[np.newaxis,np.newaxis,:]*psin[:,:,np.newaxis]
1754        B = CB[np.newaxis,np.newaxis,:]*pcos[:,:,np.newaxis]
1755        return A+B   
1756
1757    def posSawtooth(tau,Toff,slopes):
1758        Tau = (tau-Toff)%1.
1759        A = slopes[:,np.newaxis]*Tau
1760        return A
1761   
1762    def posZigZag(tau,Toff,slopes):
1763        Tau = (tau-Toff)%1.
1764        A = np.where(Tau <= 0.5,slopes[:,np.newaxis]*Tau,slopes[:,np.newaxis]*(1.-Tau))
1765        return A
1766       
1767    def DoFrac():
1768        delt2 = np.eye(2)*0.001
1769        FSC = np.ones(2,dtype='i')
1770        VFSC = np.ones(2)
1771        CSI = [np.zeros((2),dtype='i'),np.zeros(2)]
1772        if 'Crenel' in waveType:
1773            dF = fracCrenel(tau,delt2[:1],delt2[1:]).squeeze()
1774        else:
1775            dF = fracFourier(tau,nH,delt2[:1],delt2[1:]).squeeze()
1776        dFT = np.zeros_like(dF)
1777        dFTP = []
1778        for i in SdIndx:
1779            sop = Sop[i]
1780            ssop = SSop[i]
1781            ssopinv = nl.inv(ssop[0])
1782            mst = ssopinv[3][:3]
1783            epsinv = ssopinv[3][3]
1784            sdet = nl.det(sop[0])
1785            ssdet = nl.det(ssop[0])
1786            dtau = mst*(XYZ-sop[1])-epsinv*ssop[1][3]
1787            dT = 1.0
1788            if np.any(dtau%.5):
1789                dT = np.tan(np.pi*np.sum(dtau%.5))
1790            tauT = np.inner(mst,XYZ-sop[1])+epsinv*(tau-ssop[1][3])
1791            fsc = np.ones(2,dtype='i')
1792            if 'Crenel' in waveType:
1793                dFT = fracCrenel(tauT,delt2[:1],delt2[1:]).squeeze()
1794                fsc = [1,1]
1795            else:   #Fourier
1796                dFT = fracFourier(tauT,nH,delt2[:1],delt2[1:]).squeeze()
1797                dFT = nl.det(sop[0])*dFT
1798                dFT = dFT[:,np.argsort(tauT)]
1799                dFT[0] *= ssdet
1800                dFT[1] *= sdet
1801                dFTP.append(dFT)
1802           
1803                if np.any(dtau%.5) and ('1/2' in SSGData['modSymb'] or '1' in SSGData['modSymb']):
1804                    fsc = [1,1]
1805                    CSI = [[[1,0],[1,0]],[[1.,0.],[1/dT,0.]]]
1806                    FSC = np.zeros(2,dtype='i')
1807                    return CSI,dF,dFTP
1808                else:
1809                    for i in range(2):
1810                        if np.allclose(dF[i,:],dFT[i,:],atol=1.e-6):
1811                            fsc[i] = 1
1812                        else:
1813                            fsc[i] = 0
1814                    FSC &= fsc
1815                    print SSMT2text(ssop).replace(' ',''),sdet,ssdet,epsinv,fsc
1816        n = -1
1817        for i,F in enumerate(FSC):
1818            if F:
1819                n += 1
1820                CSI[0][i] = n+1
1821                CSI[1][i] = 1.0
1822       
1823        return CSI,dF,dFTP
1824       
1825    def DoXYZ():
1826        delt4 = np.ones(4)*0.001
1827        delt6 = np.eye(6)*0.001
1828        if 'Fourier' in waveType:
1829            dX = posFourier(tau,nH,delt6[:3],delt6[3:]) #+np.array(XYZ)[:,np.newaxis,np.newaxis]
1830              #3x6x12 modulated position array (X,Spos,tau)& force positive
1831            CSI = [np.zeros((6,3),dtype='i'),np.zeros((6,3))]
1832        elif waveType == 'Sawtooth':
1833            dX = posSawtooth(tau,delt4[0],delt4[1:])
1834            CSI = [np.array([[1,0,0],[2,0,0],[3,0,0],[4,0,0]]),
1835                np.array([[1.0,.0,.0],[1.0,.0,.0],[1.0,.0,.0],[1.0,.0,.0]])]
1836        elif waveType == 'ZigZag':
1837            dX = posZigZag(tau,delt4[0],delt4[1:])
1838            CSI = [np.array([[1,0,0],[2,0,0],[3,0,0],[4,0,0]]),
1839                np.array([[1.0,.0,.0],[1.0,.0,.0],[1.0,.0,.0],[1.0,.0,.0]])]
1840        XSC = np.ones(6,dtype='i')
1841        dXTP = []
1842        for i in SdIndx:
1843            sop = Sop[i]
1844            ssop = SSop[i]
1845            xsc = np.ones(6,dtype='i')
1846            ssopinv = nl.inv(ssop[0])
1847            mst = ssopinv[3][:3]
1848            epsinv = ssopinv[3][3]
1849            sdet = nl.det(sop[0])
1850            ssdet = nl.det(ssop[0])
1851            dtau = mst*(XYZ-sop[1])-epsinv*ssop[1][3]
1852            dT = 1.0
1853            if np.any(dtau%.5):
1854                dT = np.tan(np.pi*np.sum(dtau%.5))
1855            tauT = np.inner(mst,XYZ-sop[1])+epsinv*(tau-ssop[1][3])
1856            if waveType == 'Fourier':
1857                dXT = posFourier(np.sort(tauT),nH,delt6[:3],delt6[3:])   #+np.array(XYZ)[:,np.newaxis,np.newaxis]
1858            elif waveType == 'Sawtooth':
1859                dXT = posSawtooth(tauT,delt4[0],delt4[1:])+np.array(XYZ)[:,np.newaxis,np.newaxis]
1860            elif waveType == 'ZigZag':
1861                dXT = posZigZag(tauT,delt4[0],delt4[1:])+np.array(XYZ)[:,np.newaxis,np.newaxis] 
1862            dXT = np.inner(sop[0],dXT.T)    # X modulations array(3x6x49) -> array(3x49x6)
1863            dXT = np.swapaxes(dXT,1,2)      # back to array(3x6x49)
1864            dXT[:,:3,:] *= (ssdet*sdet)            # modify the sin component
1865            dXTP.append(dXT)
1866#            print np.sum(dX**2,axis=-1)*1e5
1867#            print np.sum(dXT**2,axis=-1)*1e5
1868            if waveType == 'Fourier':
1869                if np.any(dtau%.5) and ('1/2' in SSGData['modSymb'] or '1' in SSGData['modSymb']):
1870                    xsc[3:6] = 0
1871                    CSI = [[[1,0,0],[2,0,0],[3,0,0], [1,0,0],[2,0,0],[3,0,0]],
1872                        [[1.,0.,0.],[1.,0.,0.],[1.,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.]]]                   
1873                    if '(x)' in siteSym:
1874                        CSI[1][3:] = [1./dT,0.,0.],[-dT,0.,0.],[-dT,0.,0.]
1875                        if 'm' in siteSym and len(SdIndx) == 1:
1876                            CSI[1][3:] = [-dT,0.,0.],[1./dT,0.,0.],[1./dT,0.,0.]
1877                    elif '(y)' in siteSym:
1878                        CSI[1][3:] = [-dT,0.,0.],[1./dT,0.,0.],[-dT,0.,0.]
1879                        if 'm' in siteSym and len(SdIndx) == 1:
1880                            CSI[1][3:] = [1./dT,0.,0.],[-dT,0.,0.],[1./dT,0.,0.]
1881                    elif '(z)' in siteSym:
1882                        CSI[1][3:] = [-dT,0.,0.],[-dT,0.,0.],[1./dT,0.,0.]
1883                        if 'm' in siteSym and len(SdIndx) == 1:
1884                            CSI[1][3:] = [1./dT,0.,0.],[1./dT,0.,0.],[-dT,0.,0.]
1885                    elif '(xy)' in siteSym:
1886                        CSI[0] = [[1,0,0],[1,0,0],[2,0,0], [1,0,0],[1,0,0],[2,0,0]]
1887                        CSI[1][3:] = [[1./dT,0.,0.],[1./dT,0.,0.],[-dT,0.,0.]]
1888                elif '(xy)' in siteSym or '(+-0)' in siteSym:
1889                    mul = 1
1890                    if '(xy)' in siteSym:
1891                        mul = -1
1892                    if np.allclose(dX[0,0,:],dXT[1,0,:]*mul):
1893                        CSI[0][3:5] = [[11,0,0],[11,0,0]]
1894                        CSI[1][3:5] = [[1.,0,0],[mul,0,0]]
1895                        xsc[3:5] = 0
1896                    if np.allclose(dX[0,3,:],dXT[0,4,:]*mul):
1897                        CSI[0][:2] = [[12,0,0],[12,0,0]]
1898                        CSI[1][:2] = [[1.,0,0],[mul,0,0]]
1899                        xsc[:2] = 0
1900                else:
1901                    for i in range(3):
1902                        if not np.allclose(dX[i,i,:],dXT[i,i,:]):
1903                            xsc[i] = 0
1904                        if not np.allclose(dX[i,i+3,:],dXT[i,i+3,:]):
1905                            xsc[i+3] = 0
1906                XSC &= xsc
1907                print SSMT2text(ssop).replace(' ',''),sdet,ssdet,epsinv,xsc
1908        if waveType == 'Fourier':
1909            n = -1
1910            print XSC
1911            for i,X in enumerate(XSC):
1912                if X:
1913                    n += 1
1914                    CSI[0][i][0] = n+1
1915                    CSI[1][i][0] = 1.0
1916       
1917        return CSI,dX,dXTP
1918       
1919    def DoUij():
1920        tau = np.linspace(0,1,49,True)
1921        delt12 = np.eye(12)*0.0001
1922        dU = posFourier(tau,nH,delt12[:6],delt12[6:])                  #Uij modulations - 6x12x12 array
1923        CSI = [np.zeros((12,3),dtype='i'),np.zeros((12,3))]
1924        USC = np.ones(12,dtype='i')
1925        dUTP = []
1926        for i in SdIndx:
1927            sop = Sop[i]
1928            ssop = SSop[i]
1929            ssopinv = nl.inv(ssop[0])
1930            mst = ssopinv[3][:3]
1931            epsinv = ssopinv[3][3]
1932            sdet = nl.det(sop[0])
1933            ssdet = nl.det(ssop[0])
1934            dtau = mst*(XYZ-sop[1])-epsinv*ssop[1][3]
1935            dT = 1.0
1936            if np.any(dtau%.5):
1937                dT = np.tan(np.pi*np.sum(dtau%.5))
1938            tauT = np.inner(mst,XYZ-sop[1])+epsinv*(tau-ssop[1][3])
1939            usc = np.ones(12,dtype='i')
1940            dUT = posFourier(tauT,nH,delt12[:6],delt12[6:])                  #Uij modulations - 6x12x49 array
1941            dUijT = np.rollaxis(np.rollaxis(np.array(Uij2U(dUT)),3),3)    #convert dUT to 12x49x3x3
1942            dUijT = np.rollaxis(np.inner(np.inner(sop[0],dUijT),sop[0].T),3) #transform by sop - 3x3x12x49
1943            dUT = np.array(U2Uij(dUijT))    #convert to 6x12x49
1944            dUT = dUT[:,:,np.argsort(tauT)]
1945            dUTP.append(dUT)
1946            if np.any(dtau%.5) and ('1/2' in SSGData['modSymb'] or '1' in SSGData['modSymb']):
1947                CSI = [[[1,0,0],[2,0,0],[3,0,0],[4,0,0],[5,0,0],[6,0,0], 
1948                [1,0,0],[2,0,0],[3,0,0],[4,0,0],[5,0,0],[6,0,0]],
1949                [[1.,0.,0.],[1.,0.,0.],[1.,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.],
1950                [1./dT,0.,0.],[1./dT,0.,0.],[1./dT,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.]]]
1951                if 'mm2(x)' in siteSym:
1952                    CSI[1][9:] = [0.,0.,0.],[-dT,0.,0.],[0.,0.,0.]
1953                    USC = [1,1,1,0,1,0,1,1,1,0,1,0]
1954                elif '(xy)' in siteSym:
1955                    CSI[0] = [[1,0,0],[1,0,0],[2,0,0],[3,0,0],[4,0,0],[4,0,0],
1956                        [1,0,0],[1,0,0],[2,0,0],[3,0,0],[4,0,0],[4,0,0]]
1957                    CSI[1][9:] = [[1./dT,0.,0.],[-dT,0.,0.],[-dT,0.,0.]]
1958                    USC = [1,1,1,1,1,1,1,1,1,1,1,1]                             
1959                elif '(x)' in siteSym:
1960                    CSI[1][9:] = [-dT,0.,0.],[-dT,0.,0.],[1./dT,0.,0.]
1961                elif '(y)' in siteSym:
1962                    CSI[1][9:] = [-dT,0.,0.],[1./dT,0.,0.],[-dT,0.,0.]
1963                elif '(z)' in siteSym:
1964                    CSI[1][9:] = [1./dT,0.,0.],[-dT,0.,0.],[-dT,0.,0.]
1965                for i in range(6):
1966                    if not USC[i]:
1967                        CSI[0][i] = [0,0,0]
1968                        CSI[1][i] = [0.,0.,0.]
1969                        CSI[0][i+6] = [0,0,0]
1970                        CSI[1][i+6] = [0.,0.,0.]
1971            else:                       
1972                for i in range(6):
1973                    if not np.allclose(dU[i,i,:],dUT[i,i,:]*sdet):  #sin part
1974                        usc[i] = 0
1975                    if not np.allclose(dU[i,i+6,:],dUT[i,i+6,:]):   #cos part
1976                        usc[i+6] = 0
1977                if np.any(dUT[1,0,:]):
1978                    if '4/m' in siteSym:
1979                        CSI[0][6:8] = [[12,0,0],[12,0,0]]
1980                        if ssop[1][3]:
1981                            CSI[1][6:8] = [[1.,0.,0.],[-1.,0.,0.]]
1982                            usc[9] = 1
1983                        else:
1984                            CSI[1][6:8] = [[1.,0.,0.],[1.,0.,0.]]
1985                            usc[9] = 0
1986                    elif '4' in siteSym:
1987                        CSI[0][6:8] = [[12,0,0],[12,0,0]]
1988                        CSI[0][:2] = [[11,0,0],[11,0,0]]
1989                        if ssop[1][3]:
1990                            CSI[1][:2] = [[1.,0.,0.],[-1.,0.,0.]]
1991                            CSI[1][6:8] = [[1.,0.,0.],[-1.,0.,0.]]
1992                            usc[2] = 0
1993                            usc[8] = 0
1994                            usc[3] = 1
1995                            usc[9] = 1
1996                        else:
1997                            CSI[1][:2] = [[1.,0.,0.],[1.,0.,0.]]
1998                            CSI[1][6:8] = [[1.,0.,0.],[1.,0.,0.]]
1999                            usc[2] = 1
2000                            usc[8] = 1
2001                            usc[3] = 0               
2002                            usc[9] = 0
2003                    elif 'xy' in siteSym or '+-0' in siteSym:
2004                        if np.allclose(dU[0,0,:],dUT[0,1,:]*sdet):
2005                            CSI[0][4:6] = [[12,0,0],[12,0,0]]
2006                            CSI[0][6:8] = [[11,0,0],[11,0,0]]
2007                            CSI[1][4:6] = [[1.,0.,0.],[sdet,0.,0.]]
2008                            CSI[1][6:8] = [[1.,0.,0.],[sdet,0.,0.]]
2009                            usc[4:6] = 0
2010                            usc[6:8] = 0
2011                       
2012                print SSMT2text(ssop).replace(' ',''),sdet,ssdet,epsinv,usc
2013            USC &= usc
2014        print USC
2015        if not np.any(dtau%.5):
2016            n = -1
2017            for i,U in enumerate(USC):
2018                if U:
2019                    n += 1
2020                    CSI[0][i][0] = n+1
2021                    CSI[1][i][0] = 1.0
2022
2023        return CSI,dU,dUTP
2024       
2025    def DoUiso():
2026        print 'Uiso'
2027        delt4 = np.eye(4)*0.001
2028       
2029    print 'super space group: ',SSGData['SSpGrp']
2030    CSI = {'Sfrac':[[[1,0],[2,0]],[[1.,0.],[1.,0.]]],
2031        'Spos':[[[1,0,0],[2,0,0],[3,0,0], [4,0,0],[5,0,0],[6,0,0]],
2032            [[1.,0.,0.],[1.,0.,0.],[1.,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.]]],    #sin & cos
2033        'Sadp':[[[1,0,0],[2,0,0],[3,0,0],[4,0,0],[5,0,0],[6,0,0], 
2034            [7,0,0],[8,0,0],[9,0,0],[10,0,0],[11,0,0],[12,0,0]],
2035            [[1.,0.,0.],[1.,0.,0.],[1.,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.],
2036            [1.,0.,0.],[1.,0.,0.],[1.,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.]]],
2037        'Smag':[[[1,0,0],[2,0,0],[3,0,0], [4,0,0],[5,0,0],[6,0,0]],
2038            [[1.,0.,0.],[1.,0.,0.],[1.,0.,0.], [1.,0.,0.],[1.,0.,0.],[1.,0.,0.]]],}
2039    xyz = np.array(XYZ)%1.
2040    xyzt = np.array(XYZ+[0,])%1.
2041    SGOps = copy.deepcopy(SGData['SGOps'])
2042    siteSym = SytSym(XYZ,SGData)[0].strip()
2043    print 'siteSym: ',siteSym
2044    if siteSym == '1':   #"1" site symmetry
2045        if debug:
2046            return CSI,None,None,None,None
2047        else:
2048            return CSI
2049    elif siteSym == '-1':   #"-1" site symmetry
2050        CSI['Sfrac'][0] = [[1,0],[0,0]]
2051        CSI['Spos'][0] = [[1,0,0],[2,0,0],[3,0,0], [0,0,0],[0,0,0],[0,0,0]]
2052        CSI['Sadp'][0] = [[0,0,0],[0,0,0],[0,0,0],[0,0,0],[0,0,0],[0,0,0], 
2053        [1,0,0],[2,0,0],[3,0,0],[4,0,0],[5,0,0],[6,0,0]]
2054        if debug:
2055            return CSI,None,None,None,None
2056        else:
2057            return CSI
2058    SSGOps = copy.deepcopy(SSGData['SSGOps'])
2059    #expand ops to include inversions if any
2060    if SGData['SGInv']:
2061        for op,sop in zip(SGData['SGOps'],SSGData['SSGOps']):
2062            SGOps.append([-op[0],-op[1]%1.])
2063            SSGOps.append([-sop[0],-sop[1]%1.])
2064    #build set of sym ops around special position       
2065    SSop = []
2066    Sop = []
2067    Sdtau = []
2068    for iop,Op in enumerate(SGOps):         
2069        nxyz = (np.inner(Op[0],xyz)+Op[1])%1.
2070        if np.allclose(xyz,nxyz,1.e-4) and iop and MT2text(Op).replace(' ','') != '-X,-Y,-Z':
2071            SSop.append(SSGOps[iop])
2072            Sop.append(SGOps[iop])
2073            ssopinv = nl.inv(SSGOps[iop][0])
2074            mst = ssopinv[3][:3]
2075            epsinv = ssopinv[3][3]
2076            Sdtau.append(np.sum(mst*(XYZ-SGOps[iop][1])-epsinv*SSGOps[iop][1][3]))
2077    SdIndx = np.argsort(np.array(Sdtau))     # just to do in sensible order
2078    OpText =  [MT2text(s).replace(' ','') for s in Sop]         #debug?
2079    SSOpText = [SSMT2text(ss).replace(' ','') for ss in SSop]   #debug?
2080    print 'special pos super operators: ',SSOpText
2081    #setup displacement arrays
2082    tau = np.linspace(-1,1,49,True)
2083    #make modulation arrays - one parameter at a time
2084    #site fractions
2085    CSI['Sfrac'],dF,dFTP = DoFrac()
2086    #positions
2087    CSI['Spos'],dX,dXTP = DoXYZ()       
2088    #anisotropic thermal motion
2089    CSI['Sadp'],dU,dUTP = DoUij()
2090    CSI['Spos'][0] = orderParms(CSI['Spos'][0])
2091    CSI['Sadp'][0] = orderParms(CSI['Sadp'][0])           
2092    if debug:
2093        return CSI,tau,[dF,dFTP],[dX,dXTP],[dU,dUTP]
2094    else:
2095        return CSI
2096   
2097def MustrainNames(SGData):
2098    'Needs a doc string'
2099    laue = SGData['SGLaue']
2100    uniq = SGData['SGUniq']
2101    if laue in ['m3','m3m']:
2102        return ['S400','S220']
2103    elif laue in ['6/m','6/mmm','3m1']:
2104        return ['S400','S004','S202']
2105    elif laue in ['31m','3']:
2106        return ['S400','S004','S202','S211']
2107    elif laue in ['3R','3mR']:
2108        return ['S400','S220','S310','S211']
2109    elif laue in ['4/m','4/mmm']:
2110        return ['S400','S004','S220','S022']
2111    elif laue in ['mmm']:
2112        return ['S400','S040','S004','S220','S202','S022']
2113    elif laue in ['2/m']:
2114        SHKL = ['S400','S040','S004','S220','S202','S022']
2115        if uniq == 'a':
2116            SHKL += ['S013','S031','S211']
2117        elif uniq == 'b':
2118            SHKL += ['S301','S103','S121']
2119        elif uniq == 'c':
2120            SHKL += ['S130','S310','S112']
2121        return SHKL
2122    else:
2123        SHKL = ['S400','S040','S004','S220','S202','S022']
2124        SHKL += ['S310','S103','S031','S130','S301','S013']
2125        SHKL += ['S211','S121','S112']
2126        return SHKL
2127       
2128def HStrainVals(HSvals,SGData):
2129    laue = SGData['SGLaue']
2130    uniq = SGData['SGUniq']
2131    DIJ = np.zeros(6)
2132    if laue in ['m3','m3m']:
2133        DIJ[:3] = [HSvals[0],HSvals[0],HSvals[0]]
2134    elif laue in ['6/m','6/mmm','3m1','31m','3']:
2135        DIJ[:4] = [HSvals[0],HSvals[0],HSvals[1],HSvals[0]]
2136    elif laue in ['3R','3mR']:
2137        DIJ = [HSvals[0],HSvals[0],HSvals[0],HSvals[1],HSvals[1],HSvals[1]]
2138    elif laue in ['4/m','4/mmm']:
2139        DIJ[:3] = [HSvals[0],HSvals[0],HSvals[1]]
2140    elif laue in ['mmm']:
2141        DIJ[:3] = [HSvals[0],HSvals[1],HSvals[2]]
2142    elif laue in ['2/m']:
2143        DIJ[:3] = [HSvals[0],HSvals[1],HSvals[2]]
2144        if uniq == 'a':
2145            DIJ[5] = HSvals[3]
2146        elif uniq == 'b':
2147            DIJ[4] = HSvals[3]
2148        elif uniq == 'c':
2149            DIJ[3] = HSvals[3]
2150    else:
2151        DIJ = [HSvals[0],HSvals[1],HSvals[2],HSvals[3],HSvals[4],HSvals[5]]
2152    return DIJ
2153
2154def HStrainNames(SGData):
2155    'Needs a doc string'
2156    laue = SGData['SGLaue']
2157    uniq = SGData['SGUniq']
2158    if laue in ['m3','m3m']:
2159        return ['D11','eA']         #add cubic strain term
2160    elif laue in ['6/m','6/mmm','3m1','31m','3']:
2161        return ['D11','D33']
2162    elif laue in ['3R','3mR']:
2163        return ['D11','D12']
2164    elif laue in ['4/m','4/mmm']:
2165        return ['D11','D33']
2166    elif laue in ['mmm']:
2167        return ['D11','D22','D33']
2168    elif laue in ['2/m']:
2169        Dij = ['D11','D22','D33']
2170        if uniq == 'a':
2171            Dij += ['D23']
2172        elif uniq == 'b':
2173            Dij += ['D13']
2174        elif uniq == 'c':
2175            Dij += ['D12']
2176        return Dij
2177    else:
2178        Dij = ['D11','D22','D33','D12','D13','D23']
2179        return Dij
2180   
2181def MustrainCoeff(HKL,SGData):
2182    'Needs a doc string'
2183    #NB: order of terms is the same as returned by MustrainNames
2184    laue = SGData['SGLaue']
2185    uniq = SGData['SGUniq']
2186    h,k,l = HKL
2187    Strm = []
2188    if laue in ['m3','m3m']:
2189        Strm.append(h**4+k**4+l**4)
2190        Strm.append(3.0*((h*k)**2+(h*l)**2+(k*l)**2))
2191    elif laue in ['6/m','6/mmm','3m1']:
2192        Strm.append(h**4+k**4+2.0*k*h**3+2.0*h*k**3+3.0*(h*k)**2)
2193        Strm.append(l**4)
2194        Strm.append(3.0*((h*l)**2+(k*l)**2+h*k*l**2))
2195    elif laue in ['31m','3']:
2196        Strm.append(h**4+k**4+2.0*k*h**3+2.0*h*k**3+3.0*(h*k)**2)
2197        Strm.append(l**4)
2198        Strm.append(3.0*((h*l)**2+(k*l)**2+h*k*l**2))
2199        Strm.append(4.0*h*k*l*(h+k))
2200    elif laue in ['3R','3mR']:
2201        Strm.append(h**4+k**4+l**4)
2202        Strm.append(3.0*((h*k)**2+(h*l)**2+(k*l)**2))
2203        Strm.append(2.0*(h*l**3+l*k**3+k*h**3)+2.0*(l*h**3+k*l**3+l*k**3))
2204        Strm.append(4.0*(k*l*h**2+h*l*k**2+h*k*l**2))
2205    elif laue in ['4/m','4/mmm']:
2206        Strm.append(h**4+k**4)
2207        Strm.append(l**4)
2208        Strm.append(3.0*(h*k)**2)
2209        Strm.append(3.0*((h*l)**2+(k*l)**2))
2210    elif laue in ['mmm']:
2211        Strm.append(h**4)
2212        Strm.append(k**4)
2213        Strm.append(l**4)
2214        Strm.append(3.0*(h*k)**2)
2215        Strm.append(3.0*(h*l)**2)
2216        Strm.append(3.0*(k*l)**2)
2217    elif laue in ['2/m']:
2218        Strm.append(h**4)
2219        Strm.append(k**4)
2220        Strm.append(l**4)
2221        Strm.append(3.0*(h*k)**2)
2222        Strm.append(3.0*(h*l)**2)
2223        Strm.append(3.0*(k*l)**2)
2224        if uniq == 'a':
2225            Strm.append(2.0*k*l**3)
2226            Strm.append(2.0*l*k**3)
2227            Strm.append(4.0*k*l*h**2)
2228        elif uniq == 'b':
2229            Strm.append(2.0*l*h**3)
2230            Strm.append(2.0*h*l**3)
2231            Strm.append(4.0*h*l*k**2)
2232        elif uniq == 'c':
2233            Strm.append(2.0*h*k**3)
2234            Strm.append(2.0*k*h**3)
2235            Strm.append(4.0*h*k*l**2)
2236    else:
2237        Strm.append(h**4)
2238        Strm.append(k**4)
2239        Strm.append(l**4)
2240        Strm.append(3.0*(h*k)**2)
2241        Strm.append(3.0*(h*l)**2)
2242        Strm.append(3.0*(k*l)**2)
2243        Strm.append(2.0*k*h**3)
2244        Strm.append(2.0*h*l**3)
2245        Strm.append(2.0*l*k**3)
2246        Strm.append(2.0*h*k**3)
2247        Strm.append(2.0*l*h**3)
2248        Strm.append(2.0*k*l**3)
2249        Strm.append(4.0*k*l*h**2)
2250        Strm.append(4.0*h*l*k**2)
2251        Strm.append(4.0*k*h*l**2)
2252    return Strm
2253   
2254def Muiso2Shkl(muiso,SGData,cell):
2255    "this is to convert isotropic mustrain to generalized Shkls"
2256    import GSASIIlattice as G2lat
2257    A = G2lat.cell2AB(cell)[0]
2258   
2259    def minMus(Shkl,muiso,H,SGData,A):
2260        U = np.inner(A.T,H)
2261        S = np.array(MustrainCoeff(U,SGData))
2262        Sum = np.sqrt(np.sum(np.multiply(S,Shkl[:,np.newaxis]),axis=0))
2263        rad = np.sqrt(np.sum((Sum[:,np.newaxis]*H)**2,axis=1))
2264        return (muiso-rad)**2
2265       
2266    laue = SGData['SGLaue']
2267    PHI = np.linspace(0.,360.,60,True)
2268    PSI = np.linspace(0.,180.,60,True)
2269    X = np.outer(npsind(PHI),npsind(PSI))
2270    Y = np.outer(npcosd(PHI),npsind(PSI))
2271    Z = np.outer(np.ones(np.size(PHI)),npcosd(PSI))
2272    HKL = np.dstack((X,Y,Z))
2273    if laue in ['m3','m3m']:
2274        S0 = [1000.,1000.]
2275    elif laue in ['6/m','6/mmm','3m1']:
2276        S0 = [1000.,1000.,1000.]
2277    elif laue in ['31m','3']:
2278        S0 = [1000.,1000.,1000.,1000.]
2279    elif laue in ['3R','3mR']:
2280        S0 = [1000.,1000.,1000.,1000.]
2281    elif laue in ['4/m','4/mmm']:
2282        S0 = [1000.,1000.,1000.,1000.]
2283    elif laue in ['mmm']:
2284        S0 = [1000.,1000.,1000.,1000.,1000.,1000.]
2285    elif laue in ['2/m']:
2286        S0 = [1000.,1000.,1000.,0.,0.,0.,0.,0.,0.]
2287    else:
2288        S0 = [1000.,1000.,1000.,1000.,1000., 1000.,1000.,1000.,1000.,1000., 
2289            1000.,1000.,0.,0.,0.]
2290    S0 = np.array(S0)
2291    HKL = np.reshape(HKL,(-1,3))
2292    result = so.leastsq(minMus,S0,(np.ones(HKL.shape[0])*muiso,HKL,SGData,A))
2293    return result[0]
2294       
2295def SytSym(XYZ,SGData):
2296    '''
2297    Generates the number of equivalent positions and a site symmetry code for a specified coordinate and space group
2298
2299    :param XYZ: an array, tuple or list containing 3 elements: x, y & z
2300    :param SGData: from SpcGroup
2301    :Returns: a two element tuple:
2302
2303     * The 1st element is a code for the site symmetry (see GetKNsym)
2304     * The 2nd element is the site multiplicity
2305
2306    '''
2307    def PackRot(SGOps):
2308        IRT = []
2309        for ops in SGOps:
2310            M = ops[0]
2311            irt = 0
2312            for j in range(2,-1,-1):
2313                for k in range(2,-1,-1):
2314                    irt *= 3
2315                    irt += M[k][j]
2316            IRT.append(int(irt))
2317        return IRT
2318       
2319    SymName = ''
2320    Mult = 1
2321    Isym = 0
2322    if SGData['SGLaue'] in ['3','3m1','31m','6/m','6/mmm']:
2323        Isym = 1073741824
2324    Jdup = 0
2325    Xeqv = GenAtom(XYZ,SGData,True)
2326    IRT = PackRot(SGData['SGOps'])
2327    L = -1
2328    for ic,cen in enumerate(SGData['SGCen']):
2329        for invers in range(int(SGData['SGInv']+1)):
2330            for io,ops in enumerate(SGData['SGOps']):
2331                irtx = (1-2*invers)*IRT[io]
2332                L += 1
2333                if not Xeqv[L][1]:
2334                    Jdup += 1
2335                    jx = GetOprPtrName(str(irtx))
2336                    if jx[2] < 39:
2337                        Isym += 2**(jx[2]-1)
2338    if Isym == 1073741824: Isym = 0
2339    Mult = len(SGData['SGOps'])*len(SGData['SGCen'])*(int(SGData['SGInv'])+1)/Jdup
2340         
2341    return GetKNsym(str(Isym)),Mult
2342   
2343def ElemPosition(SGData):
2344    ''' Under development.
2345    Object here is to return a list of symmetry element types and locations suitable
2346    for say drawing them.
2347    So far I have the element type... getting all possible locations without lookup may be impossible!
2348    '''
2349    SymElements = []
2350    Inv = SGData['SGInv']
2351    Cen = SGData['SGCen']
2352    eleSym = {-3:['','-1'],-2:['',-6],-1:['2','-4'],0:['3','-3'],1:['4','m'],2:['6',''],3:['1','']}
2353    # get operators & expand if centrosymmetric
2354    Ops = SGData['SGOps']
2355    opM = np.array([op[0].T for op in Ops])
2356    opT = np.array([op[1] for op in Ops])
2357    if Inv:
2358        opM = np.concatenate((opM,-opM))
2359        opT = np.concatenate((opT,-opT))
2360    opMT = zip(opM,opT)
2361    for M,T in opMT[1:]:        #skip I
2362        Dt = int(nl.det(M))
2363        Tr = int(np.trace(M))
2364        Dt = -(Dt-1)/2
2365        Es = eleSym[Tr][Dt]
2366        if Dt:              #rotation-inversion
2367            I = np.eye(3)
2368            if Tr == 1:     #mirrors/glides
2369                if np.any(T):       #glide
2370                    M2 = np.inner(M,M)
2371                    MT = np.inner(M,T)+T
2372                    print 'glide',Es,MT
2373                    print M2
2374                else:               #mirror
2375                    print 'mirror',Es,T
2376                    print I-M
2377                X = [-1,-1,-1]
2378            elif Tr == -3:  # pure inversion
2379                X = np.inner(nl.inv(I-M),T)
2380                print 'inversion',Es,X
2381            else:           #other rotation-inversion
2382                M2 = np.inner(M,M)
2383                MT = np.inner(M,T)+T
2384                print 'rot-inv',Es,MT
2385                print M2
2386                X = [-1,-1,-1]
2387        else:               #rotations
2388            print 'rotation',Es
2389            X = [-1,-1,-1]
2390        #SymElements.append([Es,X])
2391       
2392    return #SymElements
2393   
2394def ApplyStringOps(A,SGData,X,Uij=[]):
2395    'Needs a doc string'
2396    SGOps = SGData['SGOps']
2397    SGCen = SGData['SGCen']
2398    Ax = A.split('+')
2399    Ax[0] = int(Ax[0])
2400    iC = 0
2401    if Ax[0] < 0:
2402        iC = 1
2403    Ax[0] = abs(Ax[0])
2404    nA = Ax[0]%100-1
2405    cA = Ax[0]/100
2406    Cen = SGCen[cA]
2407    M,T = SGOps[nA]
2408    if len(Ax)>1:
2409        cellA = Ax[1].split(',')
2410        cellA = np.array([int(a) for a in cellA])
2411    else:
2412        cellA = np.zeros(3)
2413    newX = Cen+(1-2*iC)*(np.inner(M,X)+T)+cellA
2414    if len(Uij):
2415        U = Uij2U(Uij)
2416        U = np.inner(M,np.inner(U,M).T)
2417        newUij = U2Uij(U)
2418        return [newX,newUij]
2419    else:
2420        return newX
2421       
2422def StringOpsProd(A,B,SGData):
2423    """
2424    Find A*B where A & B are in strings '-' + '100*c+n' + '+ijk'
2425    where '-' indicates inversion, c(>0) is the cell centering operator,
2426    n is operator number from SgOps and ijk are unit cell translations (each may be <0).
2427    Should return resultant string - C. SGData - dictionary using entries:
2428
2429       *  'SGCen': cell centering vectors [0,0,0] at least
2430       *  'SGOps': symmetry operations as [M,T] so that M*x+T = x'
2431
2432    """
2433    SGOps = SGData['SGOps']
2434    SGCen = SGData['SGCen']
2435    #1st split out the cell translation part & work on the operator parts
2436    Ax = A.split('+'); Bx = B.split('+')
2437    Ax[0] = int(Ax[0]); Bx[0] = int(Bx[0])
2438    iC = 0
2439    if Ax[0]*Bx[0] < 0:
2440        iC = 1
2441    Ax[0] = abs(Ax[0]); Bx[0] = abs(Bx[0])
2442    nA = Ax[0]%100-1;  nB = Bx[0]%100-1
2443    cA = Ax[0]/100;  cB = Bx[0]/100
2444    Cen = (SGCen[cA]+SGCen[cB])%1.0
2445    cC = np.nonzero([np.allclose(C,Cen) for C in SGCen])[0][0]
2446    Ma,Ta = SGOps[nA]; Mb,Tb = SGOps[nB]
2447    Mc = np.inner(Ma,Mb.T)
2448#    print Ma,Mb,Mc
2449    Tc = (np.add(np.inner(Mb,Ta)+1.,Tb))%1.0
2450#    print Ta,Tb,Tc
2451#    print [np.allclose(M,Mc)&np.allclose(T,Tc) for M,T in SGOps]
2452    nC = np.nonzero([np.allclose(M,Mc)&np.allclose(T,Tc) for M,T in SGOps])[0][0]
2453    #now the cell translation part
2454    if len(Ax)>1:
2455        cellA = Ax[1].split(',')
2456        cellA = [int(a) for a in cellA]
2457    else:
2458        cellA = [0,0,0]
2459    if len(Bx)>1:
2460        cellB = Bx[1].split(',')
2461        cellB = [int(b) for b in cellB]
2462    else:
2463        cellB = [0,0,0]
2464    cellC = np.add(cellA,cellB)
2465    C = str(((nC+1)+(100*cC))*(1-2*iC))+'+'+ \
2466        str(int(cellC[0]))+','+str(int(cellC[1]))+','+str(int(cellC[2]))
2467    return C
2468           
2469def U2Uij(U):
2470    #returns the UIJ vector U11,U22,U33,U12,U13,U23 from tensor U
2471    return [U[0][0],U[1][1],U[2][2],2.*U[0][1],2.*U[0][2],2.*U[1][2]]
2472   
2473def Uij2U(Uij):
2474    #returns the thermal motion tensor U from Uij as numpy array
2475    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]]])
2476
2477def StandardizeSpcName(spcgroup):
2478    '''Accept a spacegroup name where spaces may have not been used
2479    in the names according to the GSAS convention (spaces between symmetry
2480    for each axis) and return the space group name as used in GSAS
2481    '''
2482    rspc = spcgroup.replace(' ','').upper()
2483    # deal with rhombohedral and hexagonal setting designations
2484    rhomb = ''
2485    if rspc[-1:] == 'R':
2486        rspc = rspc[:-1]
2487        rhomb = ' R'
2488    elif rspc[-1:] == 'H': # hexagonal is assumed and thus can be ignored
2489        rspc = rspc[:-1]
2490    # look for a match in the spacegroup lists
2491    for i in spglist.values():
2492        for spc in i:
2493            if rspc == spc.replace(' ','').upper():
2494                return spc + rhomb
2495    # how about the post-2002 orthorhombic names?
2496    for i,spc in sgequiv_2002_orthorhombic:
2497        if rspc == i.replace(' ','').upper():
2498            return spc
2499    # not found
2500    return ''
2501
2502   
2503spglist = {}
2504'''A dictionary of space groups as ordered and named in the pre-2002 International
2505Tables Volume A, except that spaces are used following the GSAS convention to
2506separate the different crystallographic directions.
2507Note that the symmetry codes here will recognize many non-standard space group
2508symbols with different settings. They are ordered by Laue group
2509'''
2510spglist = {
2511    'P1' : ('P 1','P -1',), # 1-2
2512    'P2/m': ('P 2','P 21','P m','P a','P c','P n',
2513        '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
2514    'C2/m':('C 2','C m','C c','C n',
2515        'C 2/m','C 2/c','C 2/n',),
2516    'Pmmm':('P 2 2 2',
2517        'P 2 2 21','P 21 2 2','P 2 21 2',
2518        'P 21 21 2','P 2 21 21','P 21 2 21',
2519        'P 21 21 21',
2520        'P m m 2','P 2 m m','P m 2 m',
2521        'P m c 21','P 21 m a','P b 21 m','P m 21 b','P c m 21','P 21 a m',
2522        'P c c 2','P 2 a a','P b 2 b',
2523        'P m a 2','P 2 m b','P c 2 m','P m 2 a','P b m 2','P 2 c m',
2524        'P c a 21','P 21 a b','P c 21 b','P b 21 a','P b c 21','P 21 c a',
2525        'P n c 2','P 2 n a','P b 2 n','P n 2 b','P c n 2','P 2 a n',
2526        'P m n 21','P 21 m n','P n 21 m','P m 21 n','P n m 21','P 21 n m',
2527        'P b a 2','P 2 c b','P c 2 a',
2528        'P n a 21','P 21 n b','P c 21 n','P n 21 a','P b n 21','P 21 c n',
2529        'P n n 2','P 2 n n','P n 2 n',
2530        'P m m m','P n n n',
2531        'P c c m','P m a a','P b m b',
2532        'P b a n','P n c b','P c n a',
2533        'P m m a','P b m m','P m c m','P m a m','P m m b','P c m m',
2534        'P n n a','P b n n','P n c n','P n a n','P n n b','P c n n',
2535        'P m n a','P b m n','P n c m','P m a n','P n m b','P c n m',
2536        'P c c a','P b a a','P b c b','P b a b','P c c b','P c a a',
2537        'P b a m','P m c b','P c m a',
2538        'P c c n','P n a a','P b n b',
2539        'P b c m','P m c a','P b m a','P c m b','P c a m','P m a b',
2540        'P n n m','P m n n','P n m n',
2541        'P m m n','P n m m','P m n m',
2542        'P b c n','P n c a','P b n a','P c n b','P c a n','P n a b',
2543        'P b c a','P c a b',
2544        'P n m a','P b n m','P m c n','P n a m','P m n b','P c m n',
2545        ),
2546    'Cmmm':('C 2 2 21','C 2 2 2','C m m 2',
2547        'C m c 21','C c m 21','C c c 2','C m 2 m','C 2 m m',
2548        'C m 2 a','C 2 m b','C c 2 m','C 2 c m','C c 2 a','C 2 c b',
2549        'C m c m','C m c a','C c m b',
2550        'C m m m','C c c m','C m m a','C m m b','C c c a','C c c b',),
2551    'Immm':('I 2 2 2','I 21 21 21',
2552        'I m m 2','I m 2 m','I 2 m m',
2553        'I b a 2','I 2 c b','I c 2 a',
2554        'I m a 2','I 2 m b','I c 2 m','I m 2 a','I b m 2','I 2 c m',
2555        'I m m m','I b a m','I m c b','I c m a',
2556        'I b c a','I c a b',
2557        'I m m a','I b m m ','I m c m','I m a m','I m m b','I c m m',),
2558    'Fmmm':('F 2 2 2','F m m m', 'F d d d',
2559        'F m m 2','F m 2 m','F 2 m m',
2560        'F d d 2','F d 2 d','F 2 d d',),
2561    'P4/mmm':('P 4','P 41','P 42','P 43','P -4','P 4/m','P 42/m','P 4/n','P 42/n',
2562        'P 4 2 2','P 4 21 2','P 41 2 2','P 41 21 2','P 42 2 2',
2563        'P 42 21 2','P 43 2 2','P 43 21 2','P 4 m m','P 4 b m','P 42 c m',
2564        'P 42 n m','P 4 c c','P 4 n c','P 42 m c','P 42 b c','P -4 2 m',
2565        'P -4 2 c','P -4 21 m','P -4 21 c','P -4 m 2','P -4 c 2','P -4 b 2',
2566        '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',
2567        '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',
2568        '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',
2569        'P 42/n c m',),
2570    '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',
2571        'I 4 c m','I 41 m d','I 41 c d',
2572        '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',
2573        'I 41/a m d','I 41/a c d'),
2574    'R3-H':('R 3','R -3','R 3 2','R 3 m','R 3 c','R -3 m','R -3 c',),
2575    'P6/mmm': ('P 3','P 31','P 32','P -3','P 3 1 2','P 3 2 1','P 31 1 2',
2576        'P 31 2 1','P 32 1 2','P 32 2 1', 'P 3 m 1','P 3 1 m','P 3 c 1',
2577        'P 3 1 c','P -3 1 m','P -3 1 c','P -3 m 1','P -3 c 1','P 6','P 61',
2578        'P 65','P 62','P 64','P 63','P -6','P 6/m','P 63/m','P 6 2 2',
2579        'P 61 2 2','P 65 2 2','P 62 2 2','P 64 2 2','P 63 2 2','P 6 m m',
2580        'P 6 c c','P 63 c m','P 63 m c','P -6 m 2','P -6 c 2','P -6 2 m',
2581        'P -6 2 c','P 6/m m m','P 6/m c c','P 63/m c m','P 63/m m c',),
2582    'Pm3m': ('P 2 3','P 21 3','P m 3','P n 3','P a 3','P 4 3 2','P 42 3 2',
2583        'P 43 3 2','P 41 3 2','P -4 3 m','P -4 3 n','P m 3 m','P n 3 n',
2584        'P m 3 n','P n 3 m',),
2585    'Im3m':('I 2 3','I 21 3','I m -3','I a -3', 'I 4 3 2','I 41 3 2',
2586        'I -4 3 m', 'I -4 3 d','I m -3 m','I m 3 m','I a -3 d',),
2587    'Fm3m':('F 2 3','F m -3','F d -3','F 4 3 2','F 41 3 2','F -4 3 m',
2588        'F -4 3 c','F m -3 m','F m 3 m','F m -3 c','F d -3 m','F d -3 c',),
2589}
2590
2591ssdict = {}
2592'''A dictionary of superspace group symbols allowed for each entry in spglist
2593(except cubics). Monoclinics are all b-unique setting.
2594'''
2595ssdict = {
2596#1,2
2597    'P 1':['(abg)',],'P -1':['(abg)',],
2598#monoclinic - done
2599#3
2600    'P 2':['(a0g)','(a1/2g)','(0b0)','(0b0)s','(1/2b0)','(0b1/2)',],
2601#4       
2602    'P 21':['(a0g)','(0b0)','(1/2b0)','(0b1/2)',],
2603#5
2604    'C 2':['(a0g)','(0b0)','(0b0)s','(0b1/2)',],
2605#6
2606    'P m':['(a0g)','(a0g)s','(a1/2g)','(0b0)','(1/2b0)','(0b1/2)',],
2607#7
2608    'P a':['(a0g)','(a1/2g)','(0b0)','(0b1/2)',],
2609    'P c':['(a0g)','(a1/2g)','(0b0)','(1/2b0)',],
2610    'P n':['(a0g)','(a1/2g)','(0b0)','(1/2b1/2)',],
2611#8       
2612    'C m':['(a0g)','(a0g)s','(0b0)','(0b1/2)',],
2613#9       
2614    'C c':['(a0g)','(a0g)s','(0b0)',],
2615    'C n':['(a0g)','(a0g)s','(0b0)',],
2616#10       
2617    'P 2/m':['(a0g)','(a0g)0s','(a1/2g)','(0b0)','(0b0)s0','(1/2b0)','(0b1/2)',],
2618#11
2619    'P 21/m':['(a0g)','(a0g)0s','(0b0)','(0b0)s0','(1/2b0)','(0b1/2)',],
2620#12       
2621    'C 2/m':['(a0g)','(a0g)0s','(0b0)','(0b0)s0','(0b1/2)',],
2622#13
2623    'P 2/c':['(a0g)','(a0g)0s','(a1/2g)','(0b0)','(0b0)s0','(1/2b0)',],
2624    'P 2/a':['(a0g)','(a0g)0s','(a1/2g)','(0b0)','(0b0)s0','(0b1/2)',],
2625    'P 2/n':['(a0g)','(a0g)0s','(a1/2g)','(0b0)','(0b0)s0','(1/2b1/2)',],
2626#14
2627    'P 21/c':['(a0g)','(0b0)','(1/2b0)',],
2628    'P 21/a':['(a0g)','(0b0)','(0b1/2)',],
2629    'P 21/n':['(a0g)','(0b0)','(1/2b1/2)',],
2630#15
2631    'C 2/c':['(a0g)','(0b0)','(0b0)s0',],
2632    'C 2/n':['(a0g)','(0b0)','(0b0)s0',],
2633#orthorhombic
2634#16   
2635    'P 2 2 2':['(00g)','(00g)00s','(01/2g)','(1/20g)','(1/21/2g)',
2636        '(a00)','(a00)s00','(a01/2)','(a1/20)','(a1/21/2)',
2637        '(0b0)','(0b0)0s0','(1/2b0)','(0b1/2)','(1/2b1/2)',],
2638#17       
2639    'P 2 2 21':['(00g)','(01/2g)','(1/20g)','(1/21/2g)',
2640        '(a00)','(a00)s00','(a1/20)','(0b0)','(0b0)0s0','(1/2b0)',],
2641    'P 21 2 2':['(a00)','(a01/2)','(a1/20)','(a1/21/2)',
2642        '(0b0)','(0b0)0s0','(1/2b0)','(00g)','(00g)00s','(1/20g)',],
2643    'P 2 21 2':['(0b0)','(0b1/2)','(1/2b0)','(1/2b1/2)',
2644        '(00g)','(00g)00s','(1/20g)','(a00)','(a00)s00','(a1/20)',],
2645#18       
2646    'P 21 21 2':['(00g)','(00g)00s','(a00)','(a01/2)','(0b0)','(0b1/2)',],
2647    'P 2 21 21':['(a00)','(a00)s00','(0b0)','(0b1/2)','(00g)','(01/2g)',],
2648    'P 21 2 21':['(0b0)','(0b0)0s0','(00g)','(01/2g)','(a00)','(a01/2)',],
2649#19       
2650    'P 21 21 21':['(00g)','(a00)','(0b0)',],
2651#20       
2652    'C 2 2 21':['(00g)','(10g)','(01g)','(a00)','(a00)s00','(0b0)','(0b0)0s0',],
2653    'A 21 2 2':['(a00)','(a10)','(a01)','(0b0)','(0b0)0s0','(00g)','(00g)00s',],
2654    'B 2 21 2':['(0b0)','(1b0)','(0b1)','(00g)','(00g)00s','(a00)','(a00)s00',],
2655#21       
2656    'C 2 2 2':['(00g)','(00g)00s','(10g)','(10g)00s','(01g)','(01g)00s',
2657        '(a00)','(a00)s00','(a01/2)','(0b0)','(0b0)0s0','(0b1/2)',],
2658    'A 2 2 2':['(a00)','(a00)s00','(a10)','(a10)s00','(a01)','(a01)s00',
2659        '(0b0)','(0b0)0s0','(1/2b0)','(00g)','(00g)00s','(1/20g)',],
2660    'B 2 2 2':['(0b0)','(0b0)0s0','(1b0)','(1b0)0s0','(0b1)','(0b1)0s0',
2661        '(00g)','(00g)00s','(01/2g)','(a00)','(a00)s00','(a1/20)',],
2662#22       
2663    'F 2 2 2':['(00g)','(00g)00s','(10g)','(01g)',
2664        '(a00)','(a00)s00','(a10)','(a01)',
2665        '(0b0)','(0b0)0s0','(1b0)','(0b1)',],
2666#23       
2667    'I 2 2 2':['(00g)','(00g)00s','(a00)','(a00)s00','(0b0)','(0b0)0s0',],
2668#24       
2669    'I 21 21 21':['(00g)','(00g)00s','(a00)','(a00)s00','(0b0)','(0b0)0s0',],
2670#25       
2671    'P m m 2':['(00g)','(00g)s0s','(00g)0ss','(00g)ss0',
2672        '(01/2g)','(01/2g)s0s','(1/20g)','(1/20g)0ss','(1/21/2g)',
2673        '(a00)','(a00)0s0','(a1/20)','(a01/2)','(a01/2)0s0','(a1/21/2)',
2674        '(0b0)','(0b0)s00','(0b1/2)','(0b1/2)s00','(1/2b0)','(1/2b1/2)',],       
2675    'P 2 m m':['(a00)','(a00)ss0','(a00)s0s','(a00)0ss',
2676        '(a01/2)','(a01/2)ss0','(a1/20)','(a1/20)s0s','(a1/21/2)',
2677        '(0b0)','(0b0)00s','(1/2b0)','(0b1/2)','(0b1/2)00s','(1/2b1/2)',
2678        '(00g)','(00g)0s0','(01/2g)','(01/2g)0s0','(1/20g)','(1/21/2g)',],
2679    'P m 2 m':['(0b0)','(0b0)ss0','(0b0)0ss','(0b0)s0s',
2680        '(0b1/2)','(0b1/2)ss0','(1/2b0)','(1/2b0)0ss','(1/2b1/2)',
2681        '(00g)','(00g)s00','(1/20g)','(01/2g)','(01/2g)s00','(1/21/2g)',
2682        '(a00)','(a00)0s0','(a01/2)','(a01/2)0s0','(a1/20)','(a1/21/2)',],       
2683#26       
2684    'P m c 21':['(00g)','(00g)s0s','(01/2g)','(01/2g)s0s','(1/20g)','(1/21/2g)',
2685        '(a00)','(a00)0s0','(a1/20)','(0b0)','(0b0)s00','(0b1/2)',],
2686    'P 21 m a':['(a00)','(a00)ss0','(a01/2)','(a01/2)ss0','(a1/20)','(a1/21/2)',
2687        '(0b0)','(0b0)00s','(1/2b0)','(00g)','(00g)0s0','(01/2g)',],
2688    'P b 21 m':['(0b0)','(0b0)ss0','(0b1/2)','(0b1/2)ss0','(1/2b0)','(1/2b1/2)',
2689        '(00g)','(00g)s00','(1/20g)','(a00)','(a00)0s0','(a01/2)',],
2690    'P m 21 b':['(a00)','(a00)ss0','(a01/2)','(a01/2)ss0','(a1/20)','(a1/21/2)',
2691        '(00g)','(00g)0s0','(01/2g)','(0b0)','(0b0)s00','(0b1/2)',],
2692    'P c m 21':['(00g)','(00g)0ss','(1/20g)','(1/20g)0ss','(01/2g)','(1/21/2g)',
2693        '(0b0)','(0b0)s00','(1/2b0)','(a00)','(a00)0s0','(a01/2)',],
2694    'P 21 a m':['(0b0)','(0b0)ss0','(0b1/2)','(0b1/2)ss0','(1/2b0)','(1/2b1/2)',
2695        '(a00)','(a00)00s','(a1/20)','(00g)','(00g)s00','(1/20g)',],
2696#27       
2697    'P c c 2':['(00g)','(00g)s0s','(00g)0ss','(01/2g)','(1/20g)','(1/21/2g)',
2698        '(a00)','(a00)0s0','(a1/20)','(0b0)','(0b0)s00','(1/2b0)',],
2699    'P 2 a a':['(a00)','(a00)ss0','(a00)s0s','(a01/2)','(a1/20)','(a1/21/2)',
2700        '(0b0)','(0b0)00s','(0b1/2)','(00g)','(00g)0s0','(01/2g)',],
2701    'P b 2 b':['(0b0)','(0b0)0ss','(0b0)ss0','(1/2b0)','(0b1/2)','(1/2b1/2)',
2702        '(00g)','(00g)s00','(1/20g)','(a00)','(a00)00s','(a01/2)',],
2703#28       
2704    'P m a 2':['(00g)','(00g)s0s','(00g)ss0','(00g)0ss','(01/2g)','(01/2g)s0s',
2705        '(0b1/2)','(0b1/2)s00','(a01/2)','(a00)','(0b0)','(0b0)0s0','(a1/20)','(a1/21/2)'],
2706    'P 2 m b':['(a00)','(a00)s0s','(a00)ss0','(a00)0ss','(a01/2)','(a01/2)s0s',
2707        '(1/20g)','(1/20g)s00','(1/2b0)','(0b0)','(00g)','(00g)0s0','(0b1/2)','(1/2b1/2)'],
2708    'P c 2 m':['(0b0)','(0b0)s0s','(0b0)ss0','(0b0)0ss','(1/2b0)','(1/2b0)s0s',
2709        '(a1/20)','(a1/20)s00','(01/2g)','(00g)','(a00)','(a00)0s0','(1/20g)','(1/21/2g)'],
2710    'P m 2 a':['(0b0)','(0b0)s0s','(0b0)ss0','(0b0)0ss','(0b1/2)','(0b1/2)s0s',
2711        '(01/2g)','(01/2g)s00','(a1/20)','(a00)','(00g)','(00g)0s0','(a01/2)','(a1/21/2)'],
2712    'P b m 2':['(00g)','(00g)s0s','(00g)ss0','(00g)0ss','(1/20g)','(1/20g)s0s',
2713        '(a01/2)','(a01/2)s00','(0b1/2)','(0b0)','(a00)','(a00)0s0','(1/2b0)','(1/2b1/2)'],
2714    'P 2 c m':['(a00)','(a00)s0s','(a00)ss0','(a00)0ss','(a1/20)','(a1/20)s0s',
2715        '(1/2b0)','(1/2b0)s00','(1/20g)','(00g)','(0b0)','(0b0)0s0','(01/2g)','(1/21/2g)'],
2716#29       
2717    'P c a 21':['(00g)','(00g)0ss','(01/2g)','(1/20g)',
2718        '(a00)','(a00)0s0','(a1/20)','(0b0)','(0b0)s00','(1/2b0)',],
2719    'P 21 a b':['(a00)','(a00)s0s','(a01/2)','(a1/20)',
2720        '(0b0)','(0b0)00s','(0b1/2)','(00g)','(00g)0s0','(01/2g)',],
2721    'P c 21 b':['(0b0)','(0b0)ss0','(1/2b0)','(0b1/2)',
2722        '(00g)','(00g)s00','(1/20g)','(a00)','(a00)00s','(a01/2)',],
2723    'P b 21 a':['(0b0)','(0b0)0ss','(0b1/2)','(1/2b0)',
2724        '(a00)','(a00)00s','(a1/20)','(00g)','(00g)s00','(1/20g)',],
2725    'P b c 21':['(00g)','(00g)s0s','(1/20g)','(01/2g)',
2726        '(0b0)','(0b0)s00','(0b1/2)','(a00)','(a00)0s0','(a1/20)',],
2727    'P 21 c a':['(a00)','(a00)ss0','(a1/20)','(a01/2)',
2728        '(00g)','(00g)0s0','(1/20g)','(0b0)','(0b0)00s','(0b1/2)',],
2729#30       
2730    'P c n 2':['(00g)','(00g)s0s','(01/2g)','(a00)','(0b0)','(0b0)s00',
2731        '(a1/20)','(1/2b1/2)q00',],
2732    'P 2 a n':['(a00)','(a00)ss0','(a01/2)','(0b0)','(00g)','(00g)0s0',
2733        '(0b1/2)','(1/21/2g)0q0',],
2734    'P n 2 b':['(0b0)','(0b0)0ss','(1/2b0)','(00g)','(a00)','(a00)00s',
2735        '(1/20g)','(a1/21/2)00q',],
2736    'P b 2 n':['(0b0)','(0b0)ss0','(0b1/2)','(a00)','(00g)','(00g)s00',
2737        '(a01/2)','(1/21/2g)0ss',],
2738    'P n c 2':['(00g)','(00g)0ss','(1/20g)','(0b0)','(a00)','(a00)0s0',
2739        '(1/2b0)','(a1/21/2)s0s',],
2740    'P 2 n a':['(a00)','(a00)s0s','(a1/20)','(00g)','(0b0)','(0b0)00s',
2741        '(01/2g)','(1/2b1/2)ss0',],
2742#31       
2743    'P m n 21':['(00g)','(00g)s0s','(01/2g)','(01/2g)s0s','(a00)','(0b0)',
2744        '(0b0)s00','(a1/20)',],
2745    'P 21 m n':['(a00)','(a00)ss0','(a01/2)','(a01/2)ss0','(0b0)','(00g)',
2746        '(00g)0s0','(0b1/2)',],
2747    'P n 21 m':['(0b0)','(0b0)0ss','(1/2b0)','(1/2b0)0ss','(00g)','(a00)',
2748        '(a00)00s','(1/20g)',],
2749    'P m 21 n':['(0b0)','(0b0)ss0','(0b1/2)','(0b1/2)ss0','(a00)','(00g)',
2750        '(00g)s00','(a01/2)',],
2751    'P n m 21':['(00g)','(00g)0ss','(1/20g)','(1/20g)0ss','(0b0)','(a00)',
2752        '(a00)0s0','(1/2b0)',],
2753    'P 21 n m':['(a00)','(a00)s0s','(a1/20)','(a1/20)s0s','(00g)','(0b0)',
2754        '(0b0)00s','(01/2g)',],
2755#32       
2756    'P b a 2':['(00g)','(00g)s0s','(00g)0ss','(00g)ss0','(1/21/2g)qq0',
2757        '(a00)','(a01/2)','(0b0)','(0b1/2)',],
2758    'P 2 c b':['(a00)','(a00)ss0','(a00)s0s','(a00)0ss','(a1/21/2)0qq',
2759        '(0b0)','(1/2b0)','(00g)','(1/20g)',],
2760    'P c 2 a':['(0b0)','(0b0)ss0','(0b0)0ss','(0b0)s0s','(1/2b1/2)q0q',
2761        '(00g)','01/2g)','(a00)','(a1/20)',],
2762#33       
2763    'P b n 21':['(00g)','(00g)s0s','(1/21/2g)qq0','(a00)','(0b0)',],
2764    'P 21 c n':['(a00)','(a00)ss0','(a1/21/2)0qq','(0b0)','(00g)',],
2765    'P n 21 a':['(0b0)','(0b0)0ss','(1/2b1/2)q0q','(00g)','(a00)',],
2766    'P c 21 n':['(0b0)','(0b0)ss0','(1/2b1/2)q0q','(a00)','(00g)',],
2767    'P n a 21':['(00g)','(00g)0ss','(1/21/2g)qq0','(0b0)','(a00)',],
2768    'P 21 n b':['(a00)','(a00)s0s','(a1/21/2)0qq','(00g)','(0b0)',],
2769#34       
2770    'P n n 2':['(00g)','(00g)s0s','(00g)0ss','(1/21/2g)qq0',
2771        '(a00)','(a1/21/2)0q0','(a1/21/2)00q','(0b0)','(1/2b1/2)q00','(1/2b1/2)00q',],
2772    'P 2 n n':['(a00)','(a00)ss0','(a00)s0s','(a1/21/2)0qq',
2773        '(0b0)','(1/2b1/2)q00','(1/2b1/2)00q','(00g)','(1/21/2g)0q0','(1/21/2g)q00',],
2774    'P n 2 n':['(0b0)','(0b0)ss0','(0b0)0ss','(1/2b1/2)q0q',
2775        '(00g)','(1/21/2g)0q0','(1/21/2g)q00','(a00)','(a1/21/2)00q','(a1/21/2)0q0',],
2776#35       
2777    'C m m 2':['(00g)','(00g)s0s','(00g)ss0','(10g)','(10g)s0s','(10g)ss0',
2778        '(0b0)','(0b0)s00','(0b1/2)','(0b1/2)s00',],
2779    'A 2 m m':['(a00)','(a00)ss0','(a00)0ss','(a10)','(a10)ss0','(a10)0ss',
2780        '(00g)','(00g)0s0','(1/20g)','(1/20g)0s0',],
2781    'B m 2 m':['(0b0)','(0b0)0ss','(0b0)s0s','(0b1)','(0b1)0ss','(0b1)s0s',
2782        '(a00)','(a00)00s','(a1/20)','(a1/20)00s',],
2783#36
2784    'C m c 21':['(00g)','(00g)s0s','(10g)','(10g)s0s','(a00)','(a00)0s0','(0b0)','(0b0)s00',],
2785    'A 21 m a':['(a00)','(a00)ss0','(a10)','(a10)ss0','(0b0)','(0b0)00s','(00g)','(00g)0s0',],
2786    'B m 21 b':['(0b0)','(0b0)ss0','(1b0)','(1b0)ss0','(a00)','(a00)00s','(00g)','(00g)s00',],
2787    'B b 21 m':['(0b0)','(0b0)0ss','(0b1)','(0b1)ss0','(a00)','(a00)00s','(00g)','(00g)s00',],
2788    'C c m 21':['(00g)','(00g)0ss','(01g)','(01g)0ss','(a00)','(a00)0s0','(0b0)','(0b0)s00',],
2789    'A 21 a m':['(a00)','(a00)s0s','(a01)','(a01)s0s','(0b0)','(0b0)00s','(00g)','(00g)0s0',],
2790#37
2791    'C c c 2':['(00g)','(00g)s0s','(00g)0ss','(10g)','(10g)s0s','(10g)0ss','(01g)','(01g)s0s','(01g)0ss',
2792        '(a00)','(a00)0s0','(0b0)','(0b0)s00',],
2793    'A 2 a a':['(a00)','(a00)ss0','(a00)s0s','(a10)','(a10)ss0','(a10)ss0','(a01)','(a01)ss0','(a01)ss0',
2794        '(0b0)','(0b0)00s','(00g)','(00g)0s0',],
2795    'B b 2 b':['(0b0)','(0b0)0ss','(0b0)ss0','(0b1)','(0b1)0ss','(0b1)ss0','(1b0)','(1b0)0ss','(1b0)ss0',
2796        '(a00)','(a00)00s','(00g)','(00g)s00',],
2797#38
2798    'A m m 2':['(a00)','(a00)0s0','(a10)','(a10)0s0','(00g)','(00g)0s0',
2799        '(00g)ss0','(00g)0ss','(1/20g)','(1/20g)0ss','(0b0)','(0b0)s00','(1/2b0)',],
2800    'B 2 m m':['(0b0)','(0b0)00s','(0b1)','(0b1)00s','(a00)','(a00)00s',
2801        '(a00)0ss','(a00)s0s','(a1/20)','(a1/20)s0s','(00g)','(00g)0s0','(01/2g)',],
2802    'C m 2 m':['(00g)','(00g)s00','(10g)','(10g)s00','(0b0)','(0b0)s00',
2803        '(0b0)s0s','(0b0)ss0','(0b1/2)','(0b1/2)ss0','(a00)','(a00)00s','(a01/2)',],
2804    'A m 2 m':['(a00)','(a00)00s','(a01)','(a01)00s','(0b0)','(0b0)00s',
2805        '(0b0)s0s','(0b0)0ss','(1/2b0)','(1/2b0)0ss','(00g)','(00g)s00','(1/20g)',],
2806    'B m m 2':['(0b0)','(0b0)s00','(0b1)','(0b1)s00','(a00)','(a00)0s0',
2807        '(a00)0ss','(a00)ss0','(01/2g)','(01/2g)s0s','(a00)','(a00)0s0','(a1/20)',],
2808    'C 2 m m':['(00g)','(00g)0s0','(10g)','(10g)0s0','(00g)','(00g)s00',
2809        '(0b0)s0s','(0b0)0ss','(a01/2)','(a01/2)ss0','(0b0)','(0b0)00s','(0b1/2)',],
2810#39
2811    'A b m 2':['(a00)','(a00)0s0','(a01)','(a01)0s0','(00g)','(00g)s0s',
2812        '(00g)ss0','(00g)0ss','(1/20g)','(1/20g)0ss','(0b0)','(0b0)s00','(1/2b0)',],
2813    'B 2 c m':['(0b0)','(0b0)00s','(1b0)','(1b0)00s','(a00)','(a00)ss0',
2814        '(a00)0ss','(a00)s0s','(a1/20)','(a1/20)s0s','(00g)','(00g)0s0','(01/2g)',],
2815    'C m 2 a':['(00g)','(00g)s00','(01g)','(01g)s00','(0b0)','(0b0)0ss',
2816        '(0b0)s0s','(0b0)ss0','(0b1/2)','(0b1/2)ss0','(a00)','(a00)00s','(a01/2)',],
2817    'A c 2 m':['(a00)','(a00)00s','(a10)','(a10)00s','(0b0)','(0b0)ss0',
2818        '(0b0)s0s','(0b0)0ss','(1/2b0)','(1/2b0)0ss','(00g)','(00g)s00','(1/20g)',],
2819    'B m a 2':['(0b0)','(0b0)s00','(0b1)','(0b1)s00','(00g)','(00g)s0s',
2820        '(00g)0ss','(00g)ss0','(01/2g)','(01/2g)ss0','(a00)','(a00)00s','(a1/20)',],
2821    'C 2 m b':['(00g)','(00g)0s0','(10g)','(10g)0s0','(a00)','(a00)0ss',
2822        '(a00)ss0','(a00)s0s','(a01/2)','(a01/2)s0s','(0b0)','(0b0)0s0','(0b1/2)',],
2823#40       
2824    'A m a 2':['(a00)','(a01)','(00g)','(00g)s0s','(00g)ss0','(00g)0ss','(0b0)','(0b0)s00',],
2825    'B 2 m b':['(0b0)','(1b0)','(a00)','(a00)ss0','(a00)0ss','(a00)s0s','(00g)','(00g)0s0',],
2826    'C c 2 m':['(00g)','(01g)','(0b0)','(0b0)0ss','(0b0)s0s','(0b0)ss0','(a00)','(a00)00s',],
2827    'A m 2 a':['(a00)','(a10)','(0b0)','(0b0)ss0','(0b0)s0s','(0b0)0ss','(00g)','(00g)s00',],
2828    'B b m 2':['(0b0)','(0b1)','(00g)','(00g)0ss','(00g)ss0','(00g)s0s','(a00)','(a00)0s0',],
2829    'C 2 c m':['(00g)','(10g)','(a00)','(a00)s0s','(a00)0ss','(a00)ss0','(0b0)','(0b0)00s',],
2830#41
2831    'A b a 2':['(a00)','(a01)','(00g)','(00g)s0s','(00g)ss0','(00g)0ss','(0b0)','(0b0)s00',],
2832    'B 2 c b':['(0b0)','(1b0)','(a00)','(a00)ss0','(a00)0ss','(a00)s0s','(00g)','(00g)0s0',],
2833    'C c 2 a':['(00g)','(01g)','(0b0)','(0b0)0ss','(0b0)s0s','(0b0)ss0','(a00)','(a00)00s',],
2834    'A c 2 a':['(a00)','(a10)','(0b0)','(0b0)ss0','(0b0)s0s','(0b0)0ss','(00g)','(00g)s00',],
2835    'B b a 2':['(0b0)','(0b1)','(00g)','(00g)0ss','(00g)ss0','(00g)s0s','(a00)','(a00)0s0',],
2836    'C 2 c b':['(00g)','(10g)','(a00)','(a00)s0s','(a00)0ss','(a00)ss0','(0b0)','(0b0)00s',],
2837       
2838#42       
2839    'F m m 2':['(00g)','(00g)s0s','(00g)0ss','(00g)ss0','(10g)','(10g)ss0','(10g)s0s',
2840        '(01g)','(01g)ss0','(01g)0ss','(a00)','(a00)0s0','(a01)','(a01)0s0',
2841        '(0b0)','(0b0)s00','(0b1)','(0b1)s00',],       
2842    'F 2 m m':['(a00)','(a00)ss0','(a00)s0s','(a00)0ss','(a10)','(a10)0ss','(a10)ss0',
2843        '(a01)','(a01)0ss','(a01)s0s','(0b0)','(0b0)00s','(1b0)','(1b0)00s',
2844        '(00g)','(00g)0s0','(10g)','(10g)0s0',],
2845    'F m 2 m':['(0b0)','(0b0)0ss','(0b0)ss0','(0b0)s0s','(0b1)','(0b1)s0s','(0b1)0ss',
2846        '(1b0)','(1b0)s0s','(1b0)ss0','(00g)','(00g)s00','(01g)','(01g)s00',
2847        '(a00)','(a00)00s','(a10)','(a10)00s',],       
2848#43       
2849    'F d d 2':['(00g)','(00g)0ss','(00g)s0s','(a00)','(0b0)',],
2850    'F 2 d d':['(a00)','(a00)s0s','(a00)ss0','(00g)','(0b0)',],       
2851    'F d 2 d':['(0b0)','(0b0)0ss','(0b0)ss0','(a00)','(00g)',],
2852#44
2853    'I m m 2':['(00g)','(00g)ss0','(00g)s0s','(00g)0ss','(a00)','(a00)0s0','(0b0)','(0b0)s00',],
2854    'I 2 m m':['(a00)','(00g)0ss','(00g)ss0','(00g)s0s','(0b0)','(0b0)00s','(00g)','(00g)0s0',],
2855    'I m 2 m':['(0b0)','(0b0)s0s','(0b0)0ss','(0b0)ss0','(00g)','(00g)s00','(a00)','(a00)00s',],
2856#45       
2857    'I b a 2':['(00g)','(00g)ss0','(00g)s0s','(00g)0ss','(a00)','(a00)0s0','(0b0)','(0b0)s00',],
2858    'I 2 c b':['(0b0)','(0b0)s0s','(0b0)0ss','(0b0)ss0','(00g)','(00g)s00','(a00)','(a00)00s',],
2859    'I c 2 a':['(a00)','(00g)0ss','(00g)ss0','(00g)s0s','(0b0)','(0b0)00s','(00g)','(00g)0s0',],
2860#46       
2861    'I m a 2':['(a00)','(00g)0ss','(00g)ss0','(00g)s0s','(0b0)','(0b0)00s','(00g)','(00g)0s0',],
2862    'I 2 m b':['(0b0)','(0b0)s0s','(0b0)0ss','(0b0)ss0','(00g)','(00g)s00','(a00)','(a00)00s',],       
2863    'I c 2 m':['(a00)','(00g)0ss','(00g)ss0','(00g)s0s','(0b0)','(0b0)00s','(00g)','(00g)0s0',],
2864    'I m 2 a':['(a00)','(00g)0ss','(00g)ss0','(00g)s0s','(0b0)','(0b0)00s','(00g)','(00g)0s0',],
2865    'I b m 2':['(a00)','(00g)0ss','(00g)ss0','(00g)s0s','(0b0)','(0b0)00s','(00g)','(00g)0s0',],
2866    'I 2 c m':['(0b0)','(0b0)s0s','(0b0)0ss','(0b0)ss0','(00g)','(00g)s00','(a00)','(a00)00s',],
2867#47       
2868    'P m m m':['(00g)','(00g)s00','(00g)0s0','(00g)ss0','(01/2g)','(01/2g)s00','(1/20g)','(1/20g)s00','(1/21/2g)',
2869        '(a00)','(a00)0s0','(a00)00s','(a00)0ss','(a01/2)','(a01/2)0s0','(a1/20)','(a1/20)00s','(a1/21/2)',
2870        '(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s','(1/2b0)','(1/2b0)00s','(0b1/2)','(0b1/2)s00','(1/2b1/2)',],
2871#48 o@i qq0,0qq,q0q ->000
2872    'P n n n':['(00g)','(00g)s00','(00g)0s0','(1/21/2g)',
2873        '(a00)','(a00)0s0','(a00)00s','(a1/21/2)',
2874        '(0b0)','(0b0)s00','(0b0)00s','(1/2b1/2)',],
2875#49       
2876    'P c c m':['(00g)','(00g)s00','(00g)0s0','(01/2g)','(1/20g)','(1/21/2g)',
2877        '(a00)','(a00)0s0','(a00)00s','(a00)0ss','(a1/20)','(a1/20)00s',
2878        '(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s','(1/2b0)','(1/2b0)00s',],       
2879    'P m a a':['(a00)','(a00)0s0','(a00)00s','(a01/2)','(a1/20)','(a1/21/2)',
2880        '(0b0)','(0b0)00s','(0b0)s00','(0b0)s0s','(0b1/2)','(0b1/2)s00',
2881        '(00g)','(00g)0s0','(00g)s00','(00g)ss0','(01/2g)','(01/2g)s00',],       
2882    'P b m b':['(0b0)','(0b0)00s','(0b0)s00','(0b1/2)','(1/2b0)','(1/2b1/2)',
2883        '(00g)','(00g)s00','(00g)0s0','(00g)ss0','(1/20g)','(1/20g)0s0',
2884        '(a00)','(a00)00s','(a00)0s0','(a00)0ss','(a01/2)','(a01/2)0s0',],
2885#50 o@i qq0,0qq,q0q ->000
2886    'P b a n':['(00g)','(00g)s00','(00g)0s0','(00g)ss0','(1/21/2g)',
2887        '(a00)','(a00)0s0','(a01/2)','(0b0)','(0b0)s00','(0b1/2)',],
2888    'P n c b':['(a00)','(a00)0s0','(a00)00s','(a00)0ss','(a1/21/2)',
2889        '(0b0)','(0b0)00s','(1/2b0)','(00g)','(00g)0s0','(1/20g)',],
2890    'P c n a':['(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s','(1/2b1/2)',
2891        '(00g)','(00g)s00','(01/2g)','(a00)','(a00)00s','(a1/20)',],
2892#51       
2893    'P m m a':['(00g)','(00g)s00','(00g)ss0','(00g)0s0','(0b0)','(0b0)s00',
2894        '(0b0)s0s','(0b0)00s','(a00)','(a00)0s0','(01/2g)','(01/2g)s00',
2895        '(0b1/2)','(0b1/2)s00','(a01/2)','(a01/2)0s0','(1/2b0)','(1/2b1/2)',],
2896    'P b m m':['(a00)','(a00)0s0','(a00)0ss','(a00)00s','(00g)','(00g)0s0',
2897        '(00g)ss0','(00g)s00','(0b0)','(0b0)00s','(a01/2)','(a01/2)0s0',
2898        '(1/20g)','(1/20g)0s0','(1/2b0)','(1/2b0)00s','(01/2g)','(1/21/2g)',],
2899    'P m c m':['(0b0)','(0b0)00s','(0b0)s0s','(0b0)s00','(a00)','(a00)00s',
2900        '(a00)0ss','(a00)0s0','(00g)','(00g)s00','(1/2b0)','(1/2b0)00s',
2901        '(a1/20)','(a1/20)00s','(01/2g)','(01/2g)s00','(a01/2)','(a1/21/2)',],
2902    'P m a m':['(0b0)','(0b0)s00','(0b0)s0s','(0b0)00s','(00g)','(00g)s00',
2903        '(00g)ss0','(00g)0s0','(a00)','(a00)00s','(0b1/2)','(0b1/2)s00',
2904        '(01/2g)','(01/2g)s00','(a1/20)','(a1/20)00s','(1/20g)','(1/21/2g)',],
2905    'P m m b':['(00g)','(00g)0s0','(00g)ss0','(00g)s00','(a00)','(a00)0s0',
2906        '(a00)0ss','(a00)00s','(0b0)','(0b0)s00','(a00)','(a00)0s0',
2907        '(a01/2)','(a01/2)0s0','(0b1/2)','(0b1/2)s00','(a1/20)','(a1/21/2)',],
2908    'P c m m':['(a00)','(a00)00s','(a00)0ss','(a00)0s0','(0b0)','(0b0)00s',
2909        '(0b0)s0s','(0b0)s00','(00g)','(00g)0s0','(0b0)','(0b0)00s',
2910        '(1/2b0)','(1/2b0)00s','(1/20g)','(1/20g)0s0','(0b1/2)','(1/2b1/2)',],
2911#52   o@i qq0,0qq,q0q ->000     
2912    'P n n a':['(00g)','(00g)s00','(00g)0s0','(a00)','(a00)00s',
2913        '(0b0)','(0b0)00s','(a1/21/2)','(1/2b1/2)',],
2914    'P b n n':['(a00)','(a00)0s0','(a00)00s','(0b0)','(0b0)s00',
2915        '(00g)','(00g)s00','(1/2b1/2)','(1/21/2g)',],
2916    'P n c n':['(0b0)','(0b0)s00','(0b0)00s','(00g)','(00g)0s0',
2917        '(a00)','(a00)0s0','(1/21/2g)','(a1/21/2)',],
2918    'P n a n':['(0b0)','(0b0)s00','(0b0)00s','(00g)','(00g)0s0',
2919        '(a00)','(a00)0s0','(1/21/2g)','(a1/21/2)',],
2920    'P n n b':['(00g)','(00g)s00','(00g)0s0','(a00)','(a00)00s',
2921        '(0b0)','(0b0)00s','(a1/21/2)','(1/2b1/2)',],
2922    'P c n n':['(a00)','(a00)0s0','(a00)00s','(0b0)','(0b0)s00',
2923        '(00g)','(00g)s00','(1/2b1/2)','(1/21/2g)',],
2924#53       
2925    'P m n a':['(00g)','(00g)s00','(a00)','(a00)00s','(0b0)','(0b0)00s',
2926        '(0b0)s0s','(0b0)s00','(01/2g)','(01/2g)s00','(a1/20)',],
2927    'P b m n':['(a00)','(a00)0s0','(0b0)','(0b0)s00','(00g)','(00g)s00',
2928        '(00g)ss0','(00g)0s0','(a01/2)','(a01/2)0s0','(0b1/2)',],
2929    'P n c m':['(0b0)','(0b0)00s','(00g)','(00g)0s0','(a00)','(a00)0s0',
2930        '(a00)0ss','(a00)00s','(1/2b0)','(1/2b0)00s','(1/20g)',],
2931    'P m a n':['(0b0)','(0b0)s00','(a00)','(a00)0s0','(00g)','(00g)0s0',
2932        '(00g)ss0','(00g)s00','(0b1/2)','(0b1/2)s00','(a01/2)',],
2933    'P n m b':['(00g)','(00g)0s0','(0b0)','(0b0)00s','(a00)','(a00)00s',
2934        '(a00)0ss','(a00)0s0','(1/20g)','(1/20g)0s0','(1/2b0)',],
2935    'P c n m':['(a00)','(a00)00s','(00g)','(00g)s00','(0b0)','(0b0)s00',
2936        '(0b0)s0s','(0b0)00s','(a1/20)','(a1/20)00s','(01/2g)',],
2937#54       
2938    'P c c a':['(00g)','(00g)s00','(0b0)','(0b0)s00','(a00)','(a00)0s0',
2939        '(a00)0ss','(a00)00s','(01/2g)','(1/2b0)',],
2940    'P b a a':['(a00)','(a00)0s0','(00g)','(00g)0s0','(0b0)','(0b0)00s',
2941        '(0b0)s0s','(0b0)s00','(a01/2)','(01/2g)',],
2942    'P b c b':['(0b0)','(0b0)00s','(a00)','(a00)00s','(00g)','(00g)s00',
2943        '(00g)ss0','(00g)0s0','(1/2b0)','(a01/2)',],
2944    'P b a b':['(0b0)','(0b0)s00','(00g)','(00g)s00','(a00)','(a00)00s',
2945        '(a00)0ss','(a00)0s0','(0b1/2)','(1/20g)',],
2946    'P c c b':['(00g)','(00g)0s0','(a00)','(a00)0s0','(0b0)','(0b0)s00',
2947        '(0b0)s0s','(0b0)00s','(1/20g)','(a1/20)',],
2948    'P c a a':['(a00)','(a00)00s','(0b0)','(0b0)00s','(00g)','(00g)0s0',
2949        '(00g)ss0','(00g)s00','(a1/20)','(0b1/2)',],
2950#55       
2951    'P b a m':['(00g)','(00g)s00','(00g)0s0','(00g)ss0',
2952        '(a00)','(a00)00s','(a01/2)','(0b0)','(0b0)00s','(0b1/2)'],
2953    'P m c b':['(a00)','(a00)0s0','(a00)00s','(a00)0ss',
2954        '(0b0)','(0b0)s00','(1/2b0)','(00g)','(00g)s00','(1/20g)'],
2955    'P c m a':['(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s',
2956        '(a00)','(a00)0s0','(a1/20)','(00g)','(00g)0s0','(01/2g)'],
2957#56       
2958    'P c c n':['(00g)','(00g)s00','(00g)0s0','(a00)','(a00)0s0',
2959        '(0b0)','(0b0)s00'],
2960    'P n a a':['(a00)','(a00)0s0','(a00)00s','(0b0)','(0b0)00s',
2961        '(00g)','(00g)0s0'],
2962    'P b n b':['(0b0)','(0b0)s00','(0b0)00s','(a00)','(a00)00s',
2963        '(00g)','(00g)s00'],
2964#57       
2965    'P c a m':['(00g)','(00g)0s0','(a00)','(a00)00s','(0b0)','(0b0)s00',
2966        '(0b0)ss0','(0b0)00s','(01/2g)','(a1/20)','(a1/20)00s',],
2967    'P m a b':['(a00)','(a00)00s','(0b0)','(0b0)s00','(00g)','(00g)0s0',
2968        '(00g)s0s','(00g)s00','(a01/2)','(0b1/2)','(0b1/2)s00',],
2969    'P c m b':['(0b0)','(0b0)s00','(00g)','(00g)0s0','(a00)','(a00)00s',
2970        '(a00)0ss','(a00)0s0','(1/2b0)','(1/20g)','(1/20g)0s0',],
2971    'P b m a':['(0b0)','(0b0)00s','(a00)','(a00)0s0','(00g)','(00g)s00',
2972        '(00g)ss0','(00g)0s0','(0b1/2)','(a01/2)','(a01/2)0s0',],
2973    'P m c a':['(a00)','(a00)0s0','(00g)','(00g)s00','(0b0)','(0b0)00s',
2974        '(0b0)s0s','(0b0)s00','(a1/20)','(01/2g)','(01/2g)s00'],
2975    'P b c m':['(00g)','(00g)s00','(0b0)','(0b0)00s','(a00)','(a00)0s0',
2976        '(a00)0ss','(a00)00s','(1/20g)','(1/2b0)','(1/2b0)00s',],
2977#58       
2978    'P n n m':['(00g)','(00g)s00','(00g)0s0','(a00)',
2979        '(a00)00s','(0b0)','(0b0)00s'],
2980    'P m n n':['(00g)','(00g)s00','(a00)','(a00)0s0',
2981        '(a00)00s','(0b0)','(0b0)s00'],
2982    'P n m n':['(00g)','(00g)0s0','(a00)','(a00)0s0',
2983        '(0b0)','(0b0)s00','(0b0)00s',],
2984#59 o@i
2985    'P m m n':['(00g)','(00g)s00','(00g)0s0','(00g)ss0','(a00)','(a00)0s0',
2986        '(a01/2)','(a01/2)0s0','(0b0)','(0b0)s00','(0b1/2)','(0b1/2)s00',],
2987    'P n m m':['(a00)','(a00)0s0','(a00)00s','(a00)0ss','(00g)','(00g)0s0',
2988        '(1/20g)','(1/20g)0s0','(0b0)','(0b0)00s','(1/2b0)','(1/2b0)00s'],
2989    'P m n m':['(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s','(00g)','(00g)s00',
2990        '(01/2g)','(01/2g)s00','(a00)','(a00)00s','(a1/20)','(a1/20)00s'],
2991#60       
2992    'P b c n':['(00g)','(00g)s00','(00g)0s0','(a00)','(a00)0s0',
2993        '(a00)00s','(0b0)','(0b0)s00','(0b0)00s'],
2994    'P n c a':['(00g)','(00g)s00','(00g)0s0','(a00)','(a00)0s0',
2995        '(a00)00s','(0b0)','(0b0)s00','(0b0)00s'],
2996    'P b n a':['(00g)','(00g)s00','(00g)0s0','(a00)','(a00)0s0',
2997        '(a00)00s','(0b0)','(0b0)s00','(0b0)00s'],
2998    'P c n b':['(00g)','(00g)s00','(00g)0s0','(a00)','(a00)0s0',
2999        '(a00)00s','(0b0)','(0b0)s00','(0b0)00s'],
3000    'P c a n':['(00g)','(00g)s00','(00g)0s0','(a00)','(a00)0s0',
3001        '(a00)00s','(0b0)','(0b0)s00','(0b0)00s'],
3002    'P n a b':['(00g)','(00g)s00','(00g)0s0','(a00)','(a00)0s0',
3003        '(a00)00s','(0b0)','(0b0)s00','(0b0)00s'],
3004#61       
3005    'P b c a':['(00g)','(00g)s00','(00g)0s0','(a00)','(a00)0s0','(a00)00s',
3006        '(0b0)','(0b0)s00','(0b0)00s'],
3007    'P c a b':['(00g)','(00g)s00','(00g)0s0','(a00)','(a00)0s0','(a00)00s',
3008        '(0b0)','(0b0)s00','(0b0)00s'],
3009#62       
3010    'P n m a':['(00g)','(00g)0s0','(a00)','(a00)0s0','(0b0)','(0b0)00s'],
3011    'P b n m':['(00g)','(00g)s00','(a00)','(a00)00s','(0b0)','(0b0)00s'],
3012    'P m c n':['(00g)','(00g)s00','(a00)','(a00)0s0','(0b0)','(0b0)s00'],
3013    'P n a m':['(00g)','(00g)0s0','(a00)','(a00)00s','(0b0)','(0b0)00s'],
3014    'P m n b':['(00g)','(00g)s00','(a00)','(a00)00s','(0b0)','(0b0)s00'],
3015    'P c m n':['(00g)','(00g)0s0','(a00)','(a00)0s0','(0b0)','(0b0)s00'],
3016#63
3017    'C m c m':['(00g)','(00g)s00','(10g)','(10g)s00','(a00)','(a00)00s','(a00)0ss','(a00)0s0','(0b0)','(0b0)00s','(0b0)s0s','(0b0)s00',],
3018    'A m m a':['(a00)','(a00)0s0','(a10)','(a10)0s0','(0b0)','(0b0)s00','(0b0)s0s','(00g)00s','(00g)','(00g)s00','(00g)ss0','(00g)0s0',],
3019    'B b m m':['(0b0)','(0b0)00s','(0b1)','(0b1)00s','(00g)','(00g)0s0','(00g)ss0','(00g)s00','(a00)','(a00)0s0','(a00)0ss','(a00)00s',],
3020    'B m m b':['(0b0)','(0b0)s00','(1b0)','(1b0)s00','(a00)','(a00)0s0','(a00)0ss','(a00)00s','(00g)','(00g)0s0','(00g)ss0','(00g)s00',],
3021    'C c m m':['(00g)','(00g)0s0','(01g)','(01g)0s0','(0b0)','(0b0)00s','(0b0)s0s','(0b0)s00','(a00)','(a00)00s','(a00)0ss','(a00)0s0',],
3022    'A m a m':['(a00)','(a00)00s','(a01)','(a01)00s','(00g)','(00g)s00','(00g)ss0','(00g)0s0','(0b0)','(0b0)s00','(0b0)s0s','(0b0)00s',],
3023#64       
3024    'C m c a':['(00g)','(00g)s00','(10g)','(10g)s00','(0b0)','(0b0)00s','(0b0)s0s','(0b0)s00','(a00)','(a00)00s','(a00)0ss','(a00)0s0',],
3025    'A b m a':['(a00)','(a00)0s0','(a10)','(a10)0s0','(00g)','(00g)s00','(00g)ss0','(00g)0s0','(0b0)','(0b0)s00','(0b0)s0s','(0b0)00s',],
3026    'B b c m':['(0b0)','(0b0)00s','(0b1)','(0b1)00s','(a00)','(a00)0s0','(a00)0ss','(a00)00s','(00g)','(00g)0s0','(00g)ss0','(00g)s00',],
3027    'B m a b':['(0b0)','(0b0)s00','(1b0)','(1b0)s00','(00g)','(00g)0s0','(00g)ss0','(00g)s00','(a00)','(a00)0s0','(a00)0ss','(a00)00s',],
3028    'C c m b':['(00g)','(00g)0s0','(01g)','(01g)0s0','(a00)','(a00)00s','(a00)0ss','(a00)0s0','(0b0)','(0b0)00s','(0b0)s0s','(0b0)s00',],
3029    'A c a m':['(a00)','(a00)00s','(a01)','(a01)00s','(0b0)','(0b0)s00','(0b0)s0s','(0b0)00s','(00g)','(00g)s00','(00g)ss0','(00g)0s0',],
3030#65       
3031    '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',],
3032    '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',],
3033    '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',],
3034#66       
3035    'C c c m':['(00g)','(00g)s00','(10g)','(10g)s00','(0b0)','(0b0)00s','(0b0)s0s','(0b0)s00',],
3036    'A m m a':['(a00)','(a00)0s0','(a10)','(a10)0s0','(00g)','(00g)s00','(00g)ss0','(00g)0s0',],
3037    'B b m b':['(0b0)','(0b0)00s','(0b1)','(0b1)00s','(a00)','(a00)0s0','(a00)0ss','(a00)00s',],
3038#67       
3039    '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',],
3040    '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',],
3041    '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',],
3042    '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',],
3043    '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',],
3044    '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',],
3045#68 o@i
3046    'C c c a':['(00g)','(00g)s00','(10g)','(01g)','(10g)s00','(01g)s00',
3047        '(a00)','(a00)s00','(a00)ss0','(a00)0s0','(0b0)','(0b0)s00','(0b0)ss0','(0b0)0s0'],
3048    'A b a a':['(a00)','(a00)s00','(a10)','(a01)','(a10)s00','(a01)s00',
3049        '(0b0)','(0b0)s00','(0b0)ss0','(0b0)0s0','(00g)','(00g)s00','(00g)ss0','(00g)0s0'],
3050    'B b c b':['(0b0)','(0b0)s00','(0b1)','(1b0)','(0b1)s00','(1b0)s00',
3051        '(00g)','(00g)s00','(00g)ss0','(0b0)0s0','(a00)','(a00)s00','(a00)ss0','(a00)0s0'],
3052    'B b a b':['(0b0)','(0b0)s00','(1b0)','(0b1)','(1b0)s00','(0b1)s00',
3053        '(a00)','(a00)s00','(a00)ss0','(a00)0s0','(00g)','(00g)s00','(00g)ss0','(00g)0s0'],
3054    'C c c b':['(00g)','(00g)ss0','(01g)','(10g)','(01g)s00','(10g)s00',
3055        '(0b0)','(0b0)s00','(0b0)ss0','(0b0)0s0','(a00)','(a00)s00','(a00)ss0','(a00)0s0'],
3056    'A c a a':['(a00)','(a00)ss0','(a01)','(a10)','(a01)s00','(a10)s00',
3057        '(00g)','(00g)s00','(00g)ss0','(00g)0s0','(0b0)','(0b0)s00','(0b0)ss0','(0b0)0s0'],
3058#69       
3059    'F m m m':['(00g)','(00g)s00','(00g)ss0','(a00)','(a00)s00',
3060        '(a00)ss0','(0b0)','(0b0)s00','(0b0)ss0',
3061        '(10g)','(10g)s00','(10g)ss0','(a10)','(a10)0s0',
3062        '(a10)00s','(a10)0ss','(0b1)','(0b1)s00','(0b1)00s','(0b1)s0s',
3063        '(01g)','(01g)s00','(01g)ss0','(a01)','(a01)0s0',
3064        '(a01)00s','(a01)0ss','(1b0)','(1b0)s00','(1b0)00s','(1b0)s0s'],
3065#70 o@i       
3066    'F d d d':['(00g)','(00g)s00','(a00)','(a00)s00','(0b0)','(0b0)s00'],       
3067#71
3068    'I m m m':['(00g)','(00g)s00','(00g)ss0','(a00)','(a00)0s0',
3069        '(a00)ss0','(0b0)','(0b0)s00','(0b0)ss0'],
3070#72       
3071    'I b a m':['(00g)','(00g)s00','(00g)0s0','(00g)ss0','(a00)','(a00)0s0',
3072        '(a00)00s','(a00)0ss','(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s'],
3073    'I m c b':['(a00)','(a00)0s0','(a00)00s','(a00)0ss','(0b0)','(0b0)00s',
3074        '(0b0)s00','(0b0)s0s','(00g)','(00g)0s0','(00g)s00','(00g)ss0'],
3075    'I c m a':['(0b0)','(0b0)00s','(0b0)s00','(0b0)s0s','(00g)','(00g)s00',
3076        '(00g)0s0','(00g)ss0','(a00)','(a00)00s','(a00)0s0','(a00)0ss'],
3077#73       
3078    'I b c a':['(00g)','(00g)s00','(00g)0s0','(00g)ss0','(a00)','(a00)0s0',
3079        '(a00)00s','(a00)0ss','(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s'],
3080    'I c a b':['(00g)','(00g)s00','(00g)0s0','(00g)ss0','(a00)','(a00)0s0',
3081        '(a00)00s','(a00)0ss','(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s'],
3082#74       
3083    'I m m a':['(00g)','(00g)s00','(00g)0s0','(00g)ss0','(a00)','(a00)0s0',
3084        '(a00)00s','(a00)0ss','(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s'],
3085    'I b m m':['(00g)','(00g)s00','(00g)0s0','(00g)ss0','(a00)','(a00)0s0',
3086        '(a00)00s','(a00)0ss','(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s'],
3087    'I m c 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 a 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 m b':['(00g)','(00g)s00','(00g)0s0','(00g)ss0','(a00)','(a00)0s0',
3092        '(a00)00s','(a00)0ss','(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s'],
3093    'I c m m':['(00g)','(00g)s00','(00g)0s0','(00g)ss0','(a00)','(a00)0s0',
3094        '(a00)00s','(a00)0ss','(0b0)','(0b0)s00','(0b0)00s','(0b0)s0s'],
3095#tetragonal - done & checked
3096#75
3097    'P 4':['(00g)','(00g)q','(00g)s','(1/21/2g)','(1/21/2g)q',],
3098#76
3099    'P 41':['(00g)','(1/21/2g)',],
3100#77
3101    'P 42':['(00g)','(00g)q','(1/21/2g)','(1/21/2g)q',],
3102#78
3103    'P 43':['(00g)','(1/21/2g)',],
3104#79
3105    'I 4':['(00g)','(00g)q','(00g)s',],
3106#80
3107    'I 41':['(00g)','(00g)q',],
3108#81
3109    'P -4':['(00g)','(1/21/2g)',],
3110#82
3111    'I -4':['(00g)',],
3112#83
3113    'P 4/m':['(00g)','(00g)s0','(1/21/2g)',],
3114#84
3115    'P 42/m':['(00g)','(1/21/2g)',],
3116#85 o@i q0 -> 00
3117    'P 4/n':['(00g)','(00g)s0','(1/21/2g)',], #q0?
3118#86 o@i q0 -> 00
3119    'P 42/n':['(00g)','(1/21/2g)',],      #q0?
3120#87
3121    'I 4/m':['(00g)','(00g)s0',],
3122#88
3123    'I 41/a':['(00g)',],
3124#89
3125    'P 4 2 2':['(00g)','(00g)q00','(00g)s00','(1/21/2g)','(1/21/2g)q00',],
3126#90
3127    'P 4 21 2':['(00g)','(00g)q00','(00g)s00',],
3128#91
3129    'P 41 2 2':['(00g)','(1/21/2g)',],
3130#92
3131    'P 41 21 2':['(00g)',],
3132#93
3133    'P 42 2 2':['(00g)','(00g)q00','(1/21/2g)','(1/21/2g)q00',],
3134#94
3135    'P 42 21 2':['(00g)','(00g)q00',],
3136#95
3137    'P 43 2 2':['(00g)','(1/21/2g)',],
3138#96
3139    'P 43 21 2':['(00g)',],
3140#97
3141    'I 4 2 2':['(00g)','(00g)q00','(00g)s00',],
3142#98
3143    'I 41 2 2':['(00g)','(00g)q00',],
3144#99
3145    'P 4 m m':['(00g)','(00g)ss0','(00g)0ss','(00g)s0s','(1/21/2g)','(1/21/2g)0ss'],
3146#100
3147    'P 4 b m':['(00g)','(00g)ss0','(00g)0ss','(00g)s0s','(1/21/2g)qq0','(1/21/2g)qqs',],
3148#101
3149    'P 42 c m':['(00g)','(00g)0ss','(1/21/2g)','(1/21/2g)0ss',],
3150#102
3151    'P 42 n m':['(00g)','(00g)0ss','(1/21/2g)qq0','(1/21/2g)qqs',],
3152#103
3153    'P 4 c c':['(00g)','(00g)ss0','(1/21/2g)',],
3154#104
3155    'P 4 n c':['(00g)','(00g)ss0','(1/21/2g)qq0',],
3156#105
3157    'P 42 m c':['(00g)','(00g)ss0','(1/21/2g)',],
3158#106
3159    'P 42 b c':['(00g)','(00g)ss0','(1/21/2g)qq0',],
3160#107
3161    'I 4 m m':['(00g)','(00g)ss0','(00g)0ss','(00g)s0s',],
3162#108
3163    'I 4 c m':['(00g)','(00g)ss0','(00g)0ss','(00g)s0s',],
3164#109
3165    'I 41 m d':['(00g)','(00g)ss0',],
3166#110
3167    'I 41 c d':['(00g)','(00g)ss0',],
3168#111
3169    'P -4 2 m':['(00g)','(00g)0ss','(1/21/2g)','(1/21/2g)0ss',],
3170#112
3171    'P -4 2 c':['(00g)','(1/21/2g)',],
3172#113
3173    'P -4 21 m':['(00g)','(00g)0ss',],
3174#114
3175    'P -4 21 c':['(00g)',],
3176#115    00s -> 0ss
3177    'P -4 m 2':['(00g)','(00g)0s0','(1/21/2g)',],
3178#116
3179    'P -4 c 2':['(00g)','(1/21/2g)',],
3180#117    00s -> 0ss
3181    'P -4 b 2':['(00g)','(00g)0s0','(1/21/2g)0q0',],
3182#118
3183    'P -4 n 2':['(00g)','(1/21/2g)0q0',],
3184#119
3185    'I -4 m 2':['(00g)','(00g)0s0',],
3186#120
3187    'I -4 c 2':['(00g)','(00g)0s0',],
3188#121    00s -> 0ss
3189    'I -4 2 m':['(00g)','(00g)0ss',],
3190#122
3191    'I -4 2 d':['(00g)',],
3192#123
3193    'P 4/m m m':['(00g)','(00g)s0s0','(00g)00ss','(00g)s00s',
3194        '(1/21/2g)','(1/21/2g)s0s0','(1/21/2g)00ss','(1/21/2g)s00s',],
3195#124
3196    'P 4/m c c':['(00g)','(00g)s0s0','(1/21/2g)',],
3197#125    o@i q0q0 -> 0000, q0qs -> 00ss
3198    'P 4/n b m':['(00g)','(00g)s0s0','(00g)00ss','(00g)s00s','(1/21/2g)','(1/21/2g)00ss',],
3199#126    o@i q0q0 -> 0000
3200    'P 4/n n c':['(00g)','(00g)s0s0','(1/21/2g)',],
3201#127
3202    'P 4/m b m':['(00g)','(00g)s0s0','(00g)00ss','(00g)s00s',],
3203#128
3204    'P 4/m n c':['(00g)','(00g)s0s0',],
3205#129
3206    'P 4/n m m':['(00g)','(00g)s0s0','(00g)00ss','(00g)s00s',],
3207#130
3208    'P 4/n c c':['(00g)','(00g)s0s0',],
3209#131
3210    'P 42/m m c':['(00g)','(00g)s0s0','(1/21/2g)',],
3211#132
3212    'P 42/m c m':['(00g)','(00g)00ss','(1/21/2g)','(1/21/2g)00ss',],
3213#133    o@i q0q0 -> 0000
3214    'P 42/n b c':['(00g)','(00g)s0s0','(1/21/2g)',],
3215#134    o@i q0q0 -> 0000, q0qs -> 00ss
3216    'P 42/n n m':['(00g)','(00g)00ss','(1/21/2g)','(1/21/2g)00ss',],
3217#135
3218    'P 42/m b c':['(00g)','(00g)s0s0',],
3219#136
3220    'P 42/m n m':['(00g)','(00g)00ss',],
3221#137
3222    'P 42/n m c':['(00g)','(00g)s0s0',],
3223#138
3224    'P 42/n c m':['(00g)','(00g)00ss',],
3225#139
3226    'I 4/m m m':['(00g)','(00g)s0s0','(00g)00ss','(00g)s00s',],
3227#140
3228    'I 4/m c m':['(00g)','(00g)s0s0','(00g)00ss','(00g)s00s',],
3229#141
3230    'I 41/a m d':['(00g)','(00g)s0s0',],
3231#142
3232    'I 41/a c d':['(00g)','(00g)s0s0',],
3233    #trigonal/rhombahedral - done & checked
3234#143
3235    'P 3':['(00g)','(00g)t','(1/31/3g)',],
3236#144
3237    'P 31':['(00g)','(1/31/3g)',],
3238#145
3239    'P 32':['(00g)','(1/31/3g)',],
3240#146
3241    'R 3':['(00g)','(00g)t',],
3242#147
3243    'P -3':['(00g)','(1/31/3g)',],
3244#148
3245    'R -3':['(00g)',],
3246#149
3247    'P 3 1 2':['(00g)','(00g)t00','(1/31/3g)',],
3248#150
3249    'P 3 2 1':['(00g)','(00g)t00',],
3250#151
3251    'P 31 1 2':['(00g)','(1/31/3g)',],
3252#152
3253    'P 31 2 1':['(00g)',],
3254#153
3255    'P 32 1 2':['(00g)','(1/31/3g)',],
3256#154
3257    'P 32 2 1':['(00g)',],
3258#155
3259    'R 3 2':['(00g)','(00g)t0',],
3260#156
3261    'P 3 m 1':['(00g)','(00g)0s0',],
3262#157
3263    'P 3 1 m':['(00g)','(00g)00s','(1/31/3g)','(1/31/3g)00s',],
3264#158
3265    'P 3 c 1':['(00g)',],
3266#159
3267    'P 3 1 c':['(00g)','(1/31/3g)',],
3268#160
3269    'R 3 m':['(00g)','(00g)0s',],
3270#161
3271    'R 3 c':['(00g)',],
3272#162
3273    'P -3 1 m':['(00g)','(00g)00s','(1/31/3g)','(1/31/3g)00s',],
3274#163
3275    'P -3 1 c':['(00g)','(1/31/3g)',],
3276#164
3277    'P -3 m 1':['(00g)','(00g)0s0',],
3278#165
3279    'P -3 c 1':['(00g)',],
3280#166       
3281    'R -3 m':['(00g)','(00g)0s',],
3282#167
3283    'R -3 c':['(00g)',],
3284    #hexagonal - done & checked
3285#168
3286    'P 6':['(00g)','(00g)h','(00g)t','(00g)s',],
3287#169
3288    'P 61':['(00g)',],
3289#170
3290    'P 65':['(00g)',],
3291#171
3292    'P 62':['(00g)','(00g)h',],
3293#172
3294    'P 64':['(00g)','(00g)h',],
3295#173
3296    'P 63':['(00g)','(00g)h',],
3297#174
3298    'P -6':['(00g)',],
3299#175
3300    'P 6/m':['(00g)','(00g)s0',],
3301#176
3302    'P 63/m':['(00g)',],
3303#177
3304    'P 6 2 2':['(00g)','(00g)h00','(00g)t00','(00g)s00',],
3305#178
3306    'P 61 2 2':['(00g)',],
3307#179
3308    'P 65 2 2':['(00g)',],
3309#180
3310    'P 62 2 2':['(00g)','(00g)h00',],
3311#181
3312    'P 64 2 2':['(00g)','(00g)h00',],
3313#182
3314    'P 63 2 2':['(00g)','(00g)h00',],
3315#183
3316    'P 6 m m':['(00g)','(00g)ss0','(00g)0ss','(00g)s0s',],
3317#184
3318    'P 6 c c':['(00g)','(00g)s0s',],
3319#185
3320    'P 63 c m':['(00g)','(00g)0ss',],
3321#186
3322    'P 63 m c':['(00g)','(00g)0ss',],
3323#187
3324    'P -6 m 2':['(00g)','(00g)0s0',],
3325#188
3326    'P -6 c 2':['(00g)',],
3327#189
3328    'P -6 2 m':['(00g)','(00g)00s',],
3329#190
3330    'P -6 2 c':['(00g)',],
3331#191
3332    'P 6/m m m':['(00g)','(00g)s0s0','(00g)00ss','(00g)s00s',],
3333#192
3334    'P 6/m c c':['(00g)','(00g)s00s',],
3335#193
3336    'P 63/m c m':['(00g)','(00g)00ss',],
3337#194
3338    'P 63/m m c':['(00g)','(00g)00ss'],
3339    }
3340
3341#'A few non-standard space groups for test use'
3342nonstandard_sglist = ('P 21 1 1','P 1 21 1','P 1 1 21','R 3 r','R 3 2 h', 
3343                      'R -3 r', 'R 3 2 r','R 3 m h', 'R 3 m r',
3344                      'R 3 c r','R -3 c r','R -3 m r',),
3345
3346#A list of orthorhombic space groups that were renamed in the 2002 Volume A,
3347# along with the pre-2002 name. The e designates a double glide-plane'''
3348sgequiv_2002_orthorhombic= (('A e m 2', 'A b m 2',),
3349                            ('A e a 2', 'A b a 2',),
3350                            ('C m c e', 'C m c a',),
3351                            ('C m m e', 'C m m a',),
3352                            ('C c c e', 'C c c a'),)
3353#Use the space groups types in this order to list the symbols in the
3354#order they are listed in the International Tables, vol. A'''
3355symtypelist = ('triclinic', 'monoclinic', 'orthorhombic', 'tetragonal', 
3356               'trigonal', 'hexagonal', 'cubic')
3357
3358# self-test materials follow. Requires files in directory testinp
3359selftestlist = []
3360'''Defines a list of self-tests'''
3361selftestquiet = True
3362def _ReportTest():
3363    'Report name and doc string of current routine when ``selftestquiet`` is False'
3364    if not selftestquiet:
3365        import inspect
3366        caller = inspect.stack()[1][3]
3367        doc = eval(caller).__doc__
3368        if doc is not None:
3369            print('testing '+__file__+' with '+caller+' ('+doc+')')
3370        else:
3371            print('testing '+__file__()+" with "+caller)
3372def test0():
3373    '''self-test #0: exercise MoveToUnitCell'''
3374    _ReportTest()
3375    msg = "MoveToUnitCell failed"
3376    assert (MoveToUnitCell([1,2,3]) == [0,0,0]).all, msg
3377    assert (MoveToUnitCell([2,-1,-2]) == [0,0,0]).all, msg
3378    assert abs(MoveToUnitCell(np.array([-.1]))[0]-0.9) < 1e-6, msg
3379    assert abs(MoveToUnitCell(np.array([.1]))[0]-0.1) < 1e-6, msg
3380selftestlist.append(test0)
3381
3382def test1():
3383    '''self-test #1: SpcGroup against previous results'''
3384    #'''self-test #1: SpcGroup and SGPrint against previous results'''
3385    _ReportTest()
3386    testdir = ospath.join(ospath.split(ospath.abspath( __file__ ))[0],'testinp')
3387    if ospath.exists(testdir):
3388        if testdir not in sys.path: sys.path.insert(0,testdir)
3389    import spctestinp
3390    def CompareSpcGroup(spc, referr, refdict, reflist): 
3391        'Compare output from GSASIIspc.SpcGroup with results from a previous run'
3392        # if an error is reported, the dictionary can be ignored
3393        msg0 = "CompareSpcGroup failed on space group %s" % spc
3394        result = SpcGroup(spc)
3395        if result[0] == referr and referr > 0: return True
3396        keys = result[1].keys()
3397        #print result[1]['SpGrp']
3398        #msg = msg0 + " in list lengths"
3399        #assert len(keys) == len(refdict.keys()), msg
3400        for key in refdict.keys():
3401            if key == 'SGOps' or  key == 'SGCen':
3402                msg = msg0 + (" in key %s length" % key)
3403                assert len(refdict[key]) == len(result[1][key]), msg
3404                for i in range(len(refdict[key])):
3405                    msg = msg0 + (" in key %s level 0" % key)
3406                    assert np.allclose(result[1][key][i][0],refdict[key][i][0]), msg
3407                    msg = msg0 + (" in key %s level 1" % key)
3408                    assert np.allclose(result[1][key][i][1],refdict[key][i][1]), msg
3409            else:
3410                msg = msg0 + (" in key %s" % key)
3411                assert result[1][key] == refdict[key], msg
3412        msg = msg0 + (" in key %s reflist" % key)
3413        #for (l1,l2) in zip(reflist, SGPrint(result[1])):
3414        #    assert l2.replace('\t','').replace(' ','') == l1.replace(' ',''), 'SGPrint ' +msg
3415        # for now disable SGPrint testing, output has changed
3416        #assert reflist == SGPrint(result[1]), 'SGPrint ' +msg
3417    for spc in spctestinp.SGdat:
3418        CompareSpcGroup(spc, 0, spctestinp.SGdat[spc], spctestinp.SGlist[spc] )
3419selftestlist.append(test1)
3420
3421def test2():
3422    '''self-test #2: SpcGroup against cctbx (sgtbx) computations'''
3423    _ReportTest()
3424    testdir = ospath.join(ospath.split(ospath.abspath( __file__ ))[0],'testinp')
3425    if ospath.exists(testdir):
3426        if testdir not in sys.path: sys.path.insert(0,testdir)
3427    import sgtbxtestinp
3428    def CompareWcctbx(spcname, cctbx_in, debug=0):
3429        'Compare output from GSASIIspc.SpcGroup with results from cctbx.sgtbx'
3430        cctbx = cctbx_in[:] # make copy so we don't delete from the original
3431        spc = (SpcGroup(spcname))[1]
3432        if debug: print spc['SpGrp']
3433        if debug: print spc['SGCen']
3434        latticetype = spcname.strip().upper()[0]
3435        # lattice type of R implies Hexagonal centering", fix the rhombohedral settings
3436        if latticetype == "R" and len(spc['SGCen']) == 1: latticetype = 'P'
3437        assert latticetype == spc['SGLatt'], "Failed: %s does not match Lattice: %s" % (spcname, spc['SGLatt'])
3438        onebar = [1]
3439        if spc['SGInv']: onebar.append(-1)
3440        for (op,off) in spc['SGOps']:
3441            for inv in onebar:
3442                for cen in spc['SGCen']:
3443                    noff = off + cen
3444                    noff = MoveToUnitCell(noff)[0]
3445                    mult = tuple((op*inv).ravel().tolist())
3446                    if debug: print "\n%s: %s + %s" % (spcname,mult,noff)
3447                    for refop in cctbx:
3448                        if debug: print refop
3449                        # check the transform
3450                        if refop[:9] != mult: continue
3451                        if debug: print "mult match"
3452                        # check the translation
3453                        reftrans = list(refop[-3:])
3454                        reftrans = MoveToUnitCell(reftrans)[0]
3455                        if all(abs(noff - reftrans) < 1.e-5):
3456                            cctbx.remove(refop)
3457                            break
3458                    else:
3459                        assert False, "failed on %s:\n\t %s + %s" % (spcname,mult,noff)
3460    for key in sgtbxtestinp.sgtbx:
3461        CompareWcctbx(key, sgtbxtestinp.sgtbx[key])
3462selftestlist.append(test2)
3463
3464def test3(): 
3465    '''self-test #3: exercise SytSym (includes GetOprPtrName, GenAtom, GetKNsym)
3466     for selected space groups against info in IT Volume A '''
3467    _ReportTest()
3468    def ExerciseSiteSym (spc, crdlist):
3469        'compare site symmetries and multiplicities for a specified space group'
3470        msg = "failed on site sym test for %s" % spc
3471        (E,S) = SpcGroup(spc)
3472        assert not E, msg
3473        for t in crdlist:
3474            symb, m = SytSym(t[0],S)
3475            if symb.strip() != t[2].strip() or m != t[1]:
3476                print spc,t[0],m,symb,t[2]
3477            assert m == t[1]
3478            #assert symb.strip() == t[2].strip()
3479
3480    ExerciseSiteSym('p 1',[
3481            ((0.13,0.22,0.31),1,'1'),
3482            ((0,0,0),1,'1'),
3483            ])
3484    ExerciseSiteSym('p -1',[
3485            ((0.13,0.22,0.31),2,'1'),
3486            ((0,0.5,0),1,'-1'),
3487            ])
3488    ExerciseSiteSym('C 2/c',[
3489            ((0.13,0.22,0.31),8,'1'),
3490            ((0.0,.31,0.25),4,'2(y)'),
3491            ((0.25,.25,0.5),4,'-1'),
3492            ((0,0.5,0),4,'-1'),
3493            ])
3494    ExerciseSiteSym('p 2 2 2',[
3495            ((0.13,0.22,0.31),4,'1'),
3496            ((0,0.5,.31),2,'2(z)'),
3497            ((0.5,.31,0.5),2,'2(y)'),
3498            ((.11,0,0),2,'2(x)'),
3499            ((0,0.5,0),1,'222'),
3500            ])
3501    ExerciseSiteSym('p 4/n',[
3502            ((0.13,0.22,0.31),8,'1'),
3503            ((0.25,0.75,.31),4,'2(z)'),
3504            ((0.5,0.5,0.5),4,'-1'),
3505            ((0,0.5,0),4,'-1'),
3506            ((0.25,0.25,.31),2,'4(001)'),
3507            ((0.25,.75,0.5),2,'-4(001)'),
3508            ((0.25,.75,0.0),2,'-4(001)'),
3509            ])
3510    ExerciseSiteSym('p 31 2 1',[
3511            ((0.13,0.22,0.31),6,'1'),
3512            ((0.13,0.0,0.833333333),3,'2(100)'),
3513            ((0.13,0.13,0.),3,'2(110)'),
3514            ])
3515    ExerciseSiteSym('R 3 c',[
3516            ((0.13,0.22,0.31),18,'1'),
3517            ((0.0,0.0,0.31),6,'3'),
3518            ])
3519    ExerciseSiteSym('R 3 c R',[
3520            ((0.13,0.22,0.31),6,'1'),
3521            ((0.31,0.31,0.31),2,'3(111)'),
3522            ])
3523    ExerciseSiteSym('P 63 m c',[
3524            ((0.13,0.22,0.31),12,'1'),
3525            ((0.11,0.22,0.31),6,'m(100)'),
3526            ((0.333333,0.6666667,0.31),2,'3m(100)'),
3527            ((0,0,0.31),2,'3m(100)'),
3528            ])
3529    ExerciseSiteSym('I a -3',[
3530            ((0.13,0.22,0.31),48,'1'),
3531            ((0.11,0,0.25),24,'2(x)'),
3532            ((0.11,0.11,0.11),16,'3(111)'),
3533            ((0,0,0),8,'-3(111)'),
3534            ])
3535selftestlist.append(test3)
3536
3537if __name__ == '__main__':
3538    # run self-tests
3539    selftestquiet = False
3540    for test in selftestlist:
3541        test()
3542    print "OK"
Note: See TracBrowser for help on using the repository browser.