source: trunk/tests/testGSASIIspc.py @ 33

Last change on this file since 33 was 33, checked in by toby, 12 years ago

fix up package paths

File size: 5.8 KB
Line 
1# Unit tests for space group codes
2import sys
3import os
4sys.path.append('..')
5import os.path as ospath
6# determine a binary path pased on the host OS and the python version, path is relative to
7# location of this file
8if sys.platform == "win32":
9    bindir = '../binwin%d.%d' % sys.version_info[0:2]
10elif sys.platform == "darwin":
11    bindir = '../binmac%d.%d' % sys.version_info[0:2]
12else:
13    bindir = '../bin'
14if ospath.exists(ospath.join(sys.path[0],bindir)) and ospath.join(sys.path[0],bindir) not in sys.path: 
15    sys.path.insert(0,ospath.join(sys.path[0],bindir))
16
17import numpy as np
18import GSASIIspc
19
20def CompareSpcGroup(spc, referr, refdict, reflist): 
21    'Compare output from GSASIIspc.SpcGroup with results from a previous run'
22    # if an error is reported, the dictionary can be ignored
23    msg = "failed on space group %s" % spc
24    result = GSASIIspc.SpcGroup(spc)
25    if result[0] == referr and referr > 0: return True
26    keys = result[1].keys()
27    #print result[1]['SpGrp']
28    assert len(keys) == len(refdict.keys()), msg
29    for key in keys:
30        #print key, type(refdict[key])
31        if key == 'SGOps':
32            assert len(refdict[key]) == len(result[1][key]), msg
33            for i in range(len(refdict[key])):
34                assert np.allclose(result[1][key][i][0],refdict[key][i][0]), msg
35                assert np.allclose(result[1][key][i][1],refdict[key][i][1]), msg
36        else:
37            assert result[1][key] == refdict[key], msg
38    assert reflist == GSASIIspc.SGPrint(result[1]), 'SGPrint ' +msg
39
40def CompareWcctbx(spcname, cctbx_in, debug=0):
41    'Compare output from GSASIIspc.SpcGroup with results from cctbx.sgtbx'
42    cctbx = cctbx_in[:] # make copy so we don't delete from the original
43    spc = (GSASIIspc.SpcGroup(spcname))[1]
44    if debug: print spc['SpGrp']
45    if debug: print spc['SGCen']
46    latticetype = spcname.strip().upper()[0]
47    # lattice type of R implies Hexagonal centering", fix the rhombohedral settings
48    if latticetype == "R" and len(spc['SGCen']) == 1: latticetype = 'P'
49    assert latticetype == spc['SGLatt'], "Failed: %s does not match Lattice: %s" % (spcname, spc['SGLatt'])
50    onebar = [1]
51    if spc['SGInv']: onebar.append(-1)
52    for (op,off) in spc['SGOps']:
53        for inv in onebar:
54            for cen in spc['SGCen']:
55                noff = off + cen
56                GSASIIspc.MoveToUnitCell(noff)
57                mult = tuple((op*inv).ravel().tolist())
58                if debug: print "\n%s: %s + %s" % (spcname,mult,noff)
59                for refop in cctbx:
60                    if debug: print refop
61                    # check the transform
62                    if refop[:9] != mult: continue
63                    if debug: print "mult match"
64                    # check the translation
65                    reftrans = list(refop[-3:])
66                    GSASIIspc.MoveToUnitCell(reftrans)
67                    if all(abs(noff - reftrans) < 1.e-5):
68                        cctbx.remove(refop)
69                        break
70                else:
71                    assert False, "failed on %s:\n\t %s + %s" % (spcname,mult,noff)
72
73def TestSiteSym (spc, crdlist):
74    'compare site symmetries and multiplicities for a specified space group'
75    msg = "failed on site sym test for %s" % spc
76    (E,S) = GSASIIspc.SpcGroup(spc)
77    assert not E, msg
78
79    for t in crdlist:
80        symb, m = GSASIIspc.SytSym(t[0],S)
81        if symb.strip() != t[2].strip() or m != t[1]:
82            print spc,t[0],m,symb
83        assert m == t[1]
84        #assert symb.strip() == t[2].strip()
85
86#import sgtbxout
87#key = 'r -3 c'
88#CompareWcctbx(key, sgtbxout.sgtbx[key],1)
89#exit()
90
91# test #0: exercise MoveToUnitCell
92import numpy as np
93msg = "MoveToUnitCell failed"
94v = [0,1,2,-1,-2]; GSASIIspc.MoveToUnitCell(v); assert v==[0,0,0,0,0], msg
95v = np.array([-.1]); GSASIIspc.MoveToUnitCell(v); assert abs(v-0.9) < 1e-6, msg
96v = np.array([.1]); GSASIIspc.MoveToUnitCell(v); assert abs(v-0.1) < 1e-6, msg
97
98# test #1: SpcGroup and SGPrint against previous results
99import SGdat
100for spc in SGdat.SGdat:
101    #print spc
102    CompareSpcGroup(spc, 0, SGdat.SGdat[spc], SGdat.SGlist[spc] )
103
104# test #2: SpcGroup against cctbx results
105import sgtbxdat
106for key in sgtbxdat.sgtbx:
107    #print key
108    CompareWcctbx(key, sgtbxdat.sgtbx[key])
109
110# test #3: exercise SytSym (includes GetOprPtrName, GenAtom, GetKNsym)
111# for selected space groups against info in IT Volume A
112TestSiteSym('p 1',[
113        ((0.13,0.22,0.31),1,'1'),
114        ((0,0,0),1,'1'),
115])
116TestSiteSym('p -1',[
117        ((0.13,0.22,0.31),2,'1'),
118        ((0,0.5,0),1,'-1'),
119])
120TestSiteSym('C 2/c',[
121        ((0.13,0.22,0.31),8,'1'),
122        ((0.0,.31,0.25),4,'2(010)'),
123        ((0.25,.25,0.5),4,'-1'),
124        ((0,0.5,0),4,'-1'),
125])
126TestSiteSym('p 2 2 2',[
127        ((0.13,0.22,0.31),4,'1'),
128        ((0,0.5,.31),2,'2(001)'),
129        ((0.5,.31,0.5),2,'2(010)'),
130        ((.11,0,0),2,'2(100)'),
131        ((0,0.5,0),1,'222'),
132])
133TestSiteSym('p 4/n',[
134        ((0.13,0.22,0.31),8,'1'),
135        ((0.25,0.75,.31),4,'2(001)'),
136        ((0.5,0.5,0.5),4,'-1'),
137        ((0,0.5,0),4,'-1'),
138        ((0.25,0.25,.31),2,'4(001)'),
139        ((0.25,.75,0.5),2,'-4(001)'),
140        ((0.25,.75,0.0),2,'-4(001)'),
141])
142TestSiteSym('p 31 2 1',[
143        ((0.13,0.22,0.31),6,'1'),
144        ((0.13,0.0,0.833333333),3,'2(100)'),
145        ((0.13,0.13,0.),3,'2(110)'),
146])
147TestSiteSym('R 3 c',[
148        ((0.13,0.22,0.31),18,'1'),
149        ((0.0,0.0,0.31),6,'3'),
150])
151TestSiteSym('R 3 c R',[
152        ((0.13,0.22,0.31),6,'1'),
153        ((0.31,0.31,0.31),2,'3(111)'),
154])
155TestSiteSym('P 63 m c',[
156        ((0.13,0.22,0.31),12,'1'),
157        ((0.11,0.22,0.31),6,'m(100)'),
158        ((0.333333,0.6666667,0.31),2,'3m(100)'),
159        ((0,0,0.31),2,'3m(100)'),
160])
161TestSiteSym('I a -3',[
162        ((0.13,0.22,0.31),48,'1'),
163        ((0.11,0,0.25),24,'2(100)'),
164        ((0.11,0.11,0.11),16,'3(111)'),
165        ((0,0,0),8,'-3(111)'),
166])
167
168print "OK"
Note: See TracBrowser for help on using the repository browser.