1# http://pyrocko.org - GPLv3
2#
3# The Pyrocko Developers, 21st Century
4# ---|P------/S----------~Lg----------
5from __future__ import absolute_import
6from pyrocko import model, response
8import logging
9import copy
10import pickle
12logger = logging.getLogger('pyrocko.io.eventdata')
15class NoRestitution(Exception):
16 pass
19class FileNotFound(Exception):
21 def __init__(self, s):
22 self.s = s
24 def __str__(self):
25 return 'File not found: %s' % self.s
28class Problems(object):
29 def __init__(self):
30 self._problems = {}
32 def add(self, kind, nslct):
33 if kind not in self._problems:
34 self._problems[kind] = set()
35 problem = self._problems[kind]
36 problem.add(nslct)
38 def dump(self, fn):
39 f = open(fn, 'wb')
40 pickle.dump(self._problems, f)
41 f.close()
43 def load(self, fn):
44 f = open(fn, 'rb')
45 self._problems = pickle.load(f)
46 f.close()
48 def mapped(self, mapping=lambda nslct: nslct[:3]):
49 p = {}
50 for kind, problem in self._problems.items():
51 nsl = set()
52 for nslct in problem:
53 nsl.add(mapping(nslct))
54 p[kind] = nsl
56 return p
59class EventDataAccess(object):
60 '''
61 Abstract base class for event data access (see rdseed.py).
62 '''
64 def __init__(self, events=None, stations=None, datapile=None):
66 self._pile = datapile
67 self._events = events
69 if stations is None:
70 self._stations = None
71 else:
72 self._stations = {}
73 for station in stations:
74 self._stations[station.nsl()] = station
76 self._problems = Problems()
78 def get_pile(self):
79 return self._pile
81 def get_pyrocko_events(self):
82 '''
83 Extract :py:class:`model.Event` instances from the volume.
84 '''
86 if not self._events:
87 self._events = self._get_events_from_file()
88 return self._events
90 def get_pyrocko_station(self, tr, relative_event=None):
91 '''
92 Get station information for a given trace.
94 :param tr: :py:class:`trace.Trace` instance
96 :returns: :py:class:`model.station.Station` objects.
97 '''
99 self._update_stations()
100 s = copy.deepcopy(self._stations[tr.nslc_id[:3]])
101 if relative_event is not None:
102 s.set_event_relative_data(relative_event)
103 return s
105 def get_pyrocko_channel(self, tr):
106 '''
107 Get channel information for a given trace.
109 :param tr: :py:class:`trace.Trace` instance
111 :returns: :py:class:`model.station.Channel` objects.
112 '''
113 sta = self.get_station(tr)
114 return sta.get_channel(tr.channel)
116 def get_pyrocko_stations(self):
117 '''
118 Exctract a list of :py:class:`model.Station` instances.
119 '''
120 return list(self._get_stations().values())
122 def _get_stations(self, relative_event=None):
123 self._update_stations()
124 stations = copy.deepcopy(self._stations)
125 if relative_event is not None:
126 for s in stations.values():
127 s.set_event_relative_data(relative_event)
129 return stations
131 def _update_stations(self):
132 if not self._stations:
133 self._stations = {}
134 for station in self._get_stations_from_file():
135 self._stations[station.nsl()] = station
136 self._insert_channel_descriptions(self._stations)
138 def _insert_channel_descriptions(self, stations):
139 pile = self.get_pile()
140 nslc_ids = pile.gather_keys(
141 lambda tr: (tr.network, tr.station, tr.location, tr.channel))
143 for nslc in nslc_ids:
144 if nslc[:3] not in stations:
145 logger.warning(
146 'No station description for trace %s.%s.%s.%s' % nslc)
147 continue
149 sta = stations[nslc[:3]]
150 try:
151 sta.add_channel(self._get_channel_description_from_file(nslc))
152 except FileNotFound:
153 logger.warning(
154 'No channel description for trace %s.%s.%s.%s' % nslc)
156 def _get_channel_description_from_file(self, nslc):
157 return model.Channel(nslc[3], None, None, 1.)
159 def iter_traces(self, group_selector=None, trace_selector=None):
161 for traces in self.get_pile().chopper_grouped(
162 gather=lambda tr: (tr.network, tr.station, tr.location),
163 group_selector=group_selector,
164 trace_selector=trace_selector):
166 yield traces
168 def problems(self):
169 return self._problems
171 def _redundant_channel_weeder(self, redundant_channel_priorities, nslcs):
173 if redundant_channel_priorities is None:
174 return []
176 # n=network,s=station,l=location,c=channel
177 # channels by station
178 by_nsl = {}
179 for nslc in nslcs:
180 nsl = nslc[:3]
181 if nsl not in by_nsl:
182 by_nsl[nsl] = []
184 by_nsl[nsl].append(nslc)
186 # figure out what items to remove
187 to_delete = []
188 for ((h1, h2), (l1, l2)) in redundant_channel_priorities:
189 for nsl, nslcs in by_nsl.items():
190 channels = [nslc[3] for nslc in nslcs]
191 if h1 in channels and \
192 h2 in channels and \
193 l1 in channels and \
194 l2 in channels:
196 to_delete.append(nslc[:3] + (l1,))
197 to_delete.append(nslc[:3] + (l2,))
199 return to_delete
201 def get_restitution(self, tr, allowed_methods):
202 if 'integration' in allowed_methods:
203 response.IntegrationResponse()
204 else:
205 raise Exception('only "integration" restitution method is allowed')