Coverage for /builds/kinetik161/ase/ase/io/jsv.py: 83.18%

1""" 

2A module for reading and writing crystal structures from JSV 

3See http://www.jcrystal.com/steffenweber/JAVA/JSV/jsv.html 

4 

5By Jesper Friis, Jan. 2012 

6""" 

7 

8 

9import re 

10 

11import numpy as np 

12 

13import ase 

14from ase.geometry import cell_to_cellpar, cellpar_to_cell 

15from ase.spacegroup import Spacegroup, crystal 

16 

17 

18def read_jsv(f): 

19 """Reads a JSV file.""" 

20 natom = nbond = npoly = 0 

21 symbols = [] 

22 labels = [] 

23 cellpar = basis = title = bonds = poly = origin = shell_numbers = None 

24 spacegroup = 1 

25 

26 headline = f.readline().strip() 

27 

28 while True: 

29 line = f.readline() 

30 if not line: 

31 break 

32 line = line.strip() 

33 m = re.match(r'^\[([^]]+)\]\s*(.*)', line) 

34 if m is None or not line: 

35 continue 

36 tag = m.groups()[0].lower() 

37 

38 if len(m.groups()) > 1: 

39 args = m.groups()[1].split() 

40 else: 

41 args = [] 

42 

43 if tag == 'cell': 

44 cellpar = [float(x) for x in args] 

45 elif tag == 'natom': 

46 natom = int(args[0]) 

47 elif tag == 'nbond': 

48 nbond = int(args[0]) 

49 # optional margin of the bondlengths 

50 elif tag == 'npoly': 

51 npoly = int(args[0]) 

52 elif tag == 'space_group': 

53 spacegroup = Spacegroup(*tuple(int(x) for x in args)) 

54 elif tag == 'title': 

55 title = m.groups()[1] 

56 elif tag == 'atoms': 

57 symbols = [] 

58 basis = np.zeros((natom, 3), dtype=float) 

59 shell_numbers = -np.ones((natom, ), dtype=int) # float? 

60 for i in range(natom): 

61 tokens = f.readline().strip().split() 

62 labels.append(tokens[0]) 

63 symbols.append(ase.data.chemical_symbols[int(tokens[1])]) 

64 basis[i] = [float(x) for x in tokens[2:5]] 

65 if len(tokens) > 5: 

66 shell_numbers[i] = float(tokens[5]) # float? 

67 elif tag == 'bonds': 

68 for i in range(nbond): 

69 f.readline() 

70 bonds = NotImplemented 

71 elif tag == 'poly': 

72 for i in range(npoly): 

73 f.readline() 

74 poly = NotImplemented 

75 elif tag == 'origin': 

76 origin = NotImplemented 

77 else: 

78 raise ValueError(f'Unknown tag: "{tag}"') 

79 

80 if headline == 'asymmetric_unit_cell': 

81 atoms = crystal(symbols=symbols, 

82 basis=basis, 

83 spacegroup=spacegroup, 

84 cellpar=cellpar, 

85 ) 

86 elif headline == 'full_unit_cell': 

87 atoms = ase.Atoms(symbols=symbols, 

88 scaled_positions=basis, 

89 cell=cellpar_to_cell(cellpar), 

90 ) 

91 atoms.info['spacegroup'] = Spacegroup(spacegroup) 

92 elif headline == 'cartesian_cell': 

93 atoms = ase.Atoms(symbols=symbols, 

94 positions=basis, 

95 cell=cellpar_to_cell(cellpar), 

96 ) 

97 atoms.info['spacegroup'] = Spacegroup(spacegroup) 

98 else: 

99 raise ValueError(f'Invalid JSV file type: "{headline}"') 

100 

101 atoms.info['title'] = title 

102 atoms.info['labels'] = labels 

103 if bonds is not None: 

104 atoms.info['bonds'] = bonds 

105 if poly is not None: 

106 atoms.info['poly'] = poly 

107 if origin is not None: 

108 atoms.info['origin'] = origin 

109 if shell_numbers is not None: 

110 atoms.info['shell_numbers'] = shell_numbers 

111 

112 return atoms 

113 

114 

115def write_jsv(fd, atoms): 

116 """Writes JSV file.""" 

117 fd.write('asymmetric_unit_cell\n') 

118 

119 fd.write('[cell]') 

120 for v in cell_to_cellpar(atoms.cell): 

121 fd.write(' %g' % v) 

122 fd.write('\n') 

123 

124 fd.write('[natom] %d\n' % len(atoms)) 

125 fd.write('[nbond] 0\n') # FIXME 

126 fd.write('[npoly] 0\n') # FIXME 

127 

128 if 'spacegroup' in atoms.info: 

129 sg = Spacegroup(atoms.info['spacegroup']) 

130 fd.write('[space_group] %d %d\n' % (sg.no, sg.setting)) 

131 else: 

132 fd.write('[space_group] 1 1\n') 

133 

134 fd.write('[title] %s\n' % atoms.info.get('title', 'untitled')) 

135 

136 fd.write('\n') 

137 fd.write('[atoms]\n') 

138 if 'labels' in atoms.info: 

139 labels = atoms.info['labels'] 

140 else: 

141 labels = ['%s%d' % (s, i + 1) for i, s in 

142 enumerate(atoms.get_chemical_symbols())] 

143 numbers = atoms.get_atomic_numbers() 

144 scaled = atoms.get_scaled_positions() 

145 for label, n, p in zip(labels, numbers, scaled): 

146 fd.write('%-4s %2d %9.6f %9.6f %9.6f\n' 

147 % (label, n, p[0], p[1], p[2])) 

148 

149 fd.write('Label AtomicNumber x y z (repeat natom times)\n') 

150 

151 fd.write('\n') 

152 fd.write('[bonds]\n') 

153 

154 fd.write('\n') 

155 fd.write('[poly]\n') 

156 

157 fd.write('\n')