1# http://pyrocko.org - GPLv3 

2# 

3# The Pyrocko Developers, 21st Century 

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

5 

6import time 

7import re 

8import logging 

9import json 

10 

11from pyrocko import model, util 

12from pyrocko.moment_tensor import MomentTensor, symmat6 

13from .base_catalog import EarthquakeCatalog, NotFound 

14 

15from pyrocko.util import urlopen 

16 

17logger = logging.getLogger('pyrocko.client.geofon') 

18 

19km = 1000. 

20 

21 

22class Geofon(EarthquakeCatalog): 

23 ''' 

24 Access the GEOFON earthquake catalog. 

25 ''' 

26 def __init__(self, get_moment_tensors=True): 

27 self.events = {} 

28 self._get_moment_tensors = get_moment_tensors 

29 

30 def flush(self): 

31 self.events = {} 

32 

33 def iter_event_names( 

34 self, 

35 time_range=None, 

36 nmax=1000, 

37 magmin=None, 

38 latmin=-90., 

39 latmax=90., 

40 lonmin=-180., 

41 lonmax=180.): 

42 

43 logger.debug('In Geofon.iter_event_names(...)') 

44 

45 dmin = time.strftime('%Y-%m-%d', time.gmtime(time_range[0])) 

46 dmax = time.strftime('%Y-%m-%d', time.gmtime(time_range[1]+24*60*60)) 

47 

48 if magmin is None: 

49 magmin = '' 

50 else: 

51 magmin = '%g' % magmin 

52 

53 ipage = 1 

54 while True: 

55 url = ('http://geofon.gfz-potsdam.de/eqinfo/list.php?' + '&'.join([ 

56 'page=%i' % ipage, 

57 'datemin=%s' % dmin, 

58 'datemax=%s' % dmax, 

59 'latmin=%g' % latmin, 

60 'latmax=%g' % latmax, 

61 'lonmin=%g' % lonmin, 

62 'lonmax=%g' % lonmax, 

63 'magmin=%s' % magmin, 

64 'fmt=geojson', 

65 'nmax=%i' % nmax])) 

66 

67 logger.debug('Opening URL: %s' % url) 

68 page = urlopen(url).read() 

69 logger.debug('Received page (%i bytes)' % len(page)) 

70 events = self._parse_events_page(page) 

71 for ev in events: 

72 if ev.moment_tensor is True: 

73 ev.moment_tensor = self.get_mt(ev) 

74 

75 if not events: 

76 break 

77 

78 for ev in events: 

79 if time_range[0] <= ev.time and ev.time <= time_range[1]: 

80 self.events[ev.name] = ev 

81 yield ev.name 

82 

83 ipage += 1 

84 

85 def get_event(self, name): 

86 logger.debug('In Geofon.get_event("%s")' % name) 

87 

88 if name not in self.events: 

89 url = 'http://geofon.gfz-potsdam.de/eqinfo/event.php' \ 

90 '?id=%s&fmt=geojson' % name 

91 logger.debug('Opening URL: %s' % url) 

92 page = urlopen(url).read() 

93 logger.debug('Received page (%i bytes)' % len(page)) 

94 

95 try: 

96 ev = self._parse_event_page(page) 

97 self.events[name] = ev 

98 

99 except NotFound: 

100 raise NotFound(url) # reraise with url 

101 

102 ev = self.events[name] 

103 

104 if ev.moment_tensor is True: 

105 ev.moment_tensor = self.get_mt(ev) 

106 

107 return ev 

108 

109 def get_mt(self, ev): 

110 syear = time.strftime('%Y', time.gmtime(ev.time)) 

111 url = 'http://geofon.gfz-potsdam.de/data/alerts/%s/%s/mt.txt' % ( 

112 syear, ev.name) 

113 

114 try: 

115 logger.debug('Opening URL: %s' % url) 

116 page = urlopen(url).read() 

117 logger.debug('Received page (%i bytes)' % len(page)) 

118 except util.HTTPError: 

119 logger.warning('No MT found for event "%s".' % ev.name) 

120 return None 

121 

122 return self._parse_mt_page(page) 

123 

124 def _parse_mt_page(self, page): 

125 d = {} 

126 for k in 'Scale', 'Mrr', 'Mtt', 'Mpp', 'Mrt', 'Mrp', 'Mtp': 

127 r = k.encode('ascii')+br'\s*=?\s*(\S+)' 

128 m = re.search(r, page) 

129 if m: 

130 s = m.group(1).replace(b'10**', b'1e') 

131 d[k.lower()] = float(s) 

132 

133 m = symmat6(*(d[x] for x in 'mrr mtt mpp mrt mrp mtp'.split())) 

134 m *= d['scale'] 

135 mt = MomentTensor(m_up_south_east=m) 

136 

137 return mt 

138 

139 def _parse_events_page(self, page, limit=None): 

140 j = json.loads(page.decode('utf-8')) 

141 events = [] 

142 for ifeature, feature in enumerate(j['features']): 

143 ev = self._json_feature_to_event(feature) 

144 events.append(ev) 

145 if limit and ifeature + 1 == limit: 

146 break 

147 

148 return events 

149 

150 def _parse_event_page(self, page): 

151 return self._parse_events_page(page, limit=1)[0] 

152 

153 def _json_feature_to_event(self, feature): 

154 name = feature['id'] 

155 lon, lat, depth = feature['geometry']['coordinates'] 

156 depth *= 1000. 

157 properties = feature['properties'] 

158 magnitude = properties['mag'] 

159 magnitude_type = properties['magType'] 

160 region = properties['place'] 

161 tevent = util.str_to_time(properties['time'].replace('T', ' ')) 

162 

163 if ((properties.get('hasMT', 'no') == 'yes') 

164 or properties['magType'] == 'Mw') and self._get_moment_tensors: 

165 

166 moment_tensor = True # flag for caller to query MT 

167 else: 

168 moment_tensor = None 

169 

170 status = properties['status'][:1] 

171 tags = [] 

172 if status in 'AMC': 

173 tags.append('geofon_status:%s' % status) 

174 

175 category = properties.get('evtype', '') 

176 if re.match(r'^[a-zA-Z0-9]+$', category): 

177 tags.append('geofon_category:%s' % category) 

178 

179 ev = model.Event( 

180 lat=float(lat), 

181 lon=float(lon), 

182 time=tevent, 

183 name=name, 

184 depth=float(depth), 

185 magnitude=float(magnitude), 

186 magnitude_type=str(magnitude_type), 

187 region=str(region), 

188 moment_tensor=moment_tensor, 

189 catalog='GEOFON', 

190 tags=tags) 

191 

192 return ev