1# http://pyrocko.org - GPLv3 

2# 

3# The Pyrocko Developers, 21st Century 

4# ---|P------/S----------~Lg---------- 

5from __future__ import division, absolute_import 

6 

7from struct import unpack 

8import os 

9import re 

10import logging 

11 

12from pyrocko import trace 

13from pyrocko.util import reuse, ensuredirs 

14from .io_common import FileLoadError, FileSaveError 

15 

16logger = logging.getLogger('pyrocko.io.mseed') 

17 

18VALID_RECORD_LENGTHS = tuple(2**exp for exp in range(8, 20)) 

19 

20 

21class CodeTooLong(FileSaveError): 

22 pass 

23 

24 

25def iload(filename, load_data=True): 

26 from pyrocko import mseed_ext 

27 

28 have_zero_rate_traces = False 

29 try: 

30 traces = [] 

31 for tr in mseed_ext.get_traces(filename, load_data): 

32 network, station, location, channel = tr[1:5] 

33 tmin = float(tr[5])/float(mseed_ext.HPTMODULUS) 

34 tmax = float(tr[6])/float(mseed_ext.HPTMODULUS) 

35 try: 

36 deltat = reuse(1.0/float(tr[7])) 

37 except ZeroDivisionError: 

38 have_zero_rate_traces = True 

39 continue 

40 

41 ydata = tr[8] 

42 

43 traces.append(trace.Trace( 

44 network, station, location, channel, tmin, tmax, 

45 deltat, ydata)) 

46 

47 for tr in traces: 

48 yield tr 

49 

50 except (OSError, mseed_ext.MSeedError) as e: 

51 raise FileLoadError(str(e)+' (file: %s)' % filename) 

52 

53 if have_zero_rate_traces: 

54 logger.warning( 

55 'Ignoring traces with sampling rate of zero in file %s ' 

56 '(maybe LOG traces)' % filename) 

57 

58 

59def as_tuple(tr): 

60 from pyrocko import mseed_ext 

61 itmin = int(round(tr.tmin*mseed_ext.HPTMODULUS)) 

62 itmax = int(round(tr.tmax*mseed_ext.HPTMODULUS)) 

63 srate = 1.0/tr.deltat 

64 return (tr.network, tr.station, tr.location, tr.channel, 

65 itmin, itmax, srate, tr.get_ydata()) 

66 

67 

68def save(traces, filename_template, additional={}, overwrite=True, 

69 steim=1, record_length=4096): 

70 from pyrocko import mseed_ext 

71 

72 assert record_length in VALID_RECORD_LENGTHS 

73 

74 fn_tr = {} 

75 for tr in traces: 

76 for code, maxlen, val in zip( 

77 ['network', 'station', 'location', 'channel'], 

78 [2, 5, 2, 3], 

79 tr.nslc_id): 

80 

81 if len(val) > maxlen: 

82 raise CodeTooLong( 

83 '%s code too long to be stored in MSeed file: %s' % 

84 (code, val)) 

85 

86 fn = tr.fill_template(filename_template, **additional) 

87 if not overwrite and os.path.exists(fn): 

88 raise FileSaveError('File exists: %s' % fn) 

89 

90 if fn not in fn_tr: 

91 fn_tr[fn] = [] 

92 

93 fn_tr[fn].append(tr) 

94 

95 for fn, traces_thisfile in fn_tr.items(): 

96 trtups = [] 

97 traces_thisfile.sort(key=lambda a: a.full_id) 

98 for tr in traces_thisfile: 

99 trtups.append(as_tuple(tr)) 

100 

101 ensuredirs(fn) 

102 try: 

103 mseed_ext.store_traces( 

104 trtups, fn, record_length=record_length, steim=steim) 

105 except mseed_ext.MSeedError as e: 

106 raise FileSaveError( 

107 str(e) + ' (while storing traces to file \'%s\')' % fn) 

108 

109 return list(fn_tr.keys()) 

110 

111 

112tcs = {} 

113 

114 

115def detect(first512): 

116 

117 if len(first512) < 256: 

118 return False 

119 

120 rec = first512 

121 

122 try: 

123 sequence_number = int(rec[:6]) 

124 except Exception: 

125 return False 

126 if sequence_number < 0: 

127 return False 

128 

129 type_code = rec[6:7] 

130 if type_code in b'DRQM': 

131 bads = [] 

132 for sex in '<>': 

133 bad = False 

134 fmt = sex + '6s1s1s5s2s3s2s10sH2h4Bl2H' 

135 vals = unpack(fmt, rec[:48]) 

136 fmt_btime = sex + 'HHBBBBH' 

137 tvals = unpack(fmt_btime, vals[7]) 

138 if tvals[1] < 1 or tvals[1] > 367 or tvals[2] > 23 or \ 

139 tvals[3] > 59 or tvals[4] > 60 or tvals[6] > 9999: 

140 bad = True 

141 

142 bads.append(bad) 

143 

144 if all(bads): 

145 return False 

146 

147 else: 

148 if type_code not in b'VAST': 

149 return False 

150 

151 continuation_code = rec[7:8] 

152 if continuation_code != b' ': 

153 return False 

154 

155 blockette_type = rec[8:8+3].decode() 

156 if not re.match(r'^\d\d\d$', blockette_type): 

157 return False 

158 

159 try: 

160 blockette_length = int(rec[11:11+4]) 

161 except Exception: 

162 return False 

163 

164 if blockette_length < 7: 

165 return False 

166 

167 return True