Browse Source

squirrel: stationxml export

release
Sebastian Heimann 3 weeks ago
parent
commit
9e5e770d68
  1. 55
      src/io/stationxml.py
  2. 185
      src/squirrel/base.py
  3. 7
      src/squirrel/cache.py
  4. 18
      src/squirrel/model.py
  5. 6
      src/squirrel/tool/commands/coverage.py
  6. 23
      src/squirrel/tool/commands/jackseis.py
  7. 60
      src/squirrel/tool/commands/stationxml.py

55
src/io/stationxml.py

@ -272,7 +272,33 @@ def tts(t):
if t is None:
return '?'
else:
return util.tts(t, format='%Y-%m-%d')
return util.tts(t, format='%Y-%m-%d %H:%M:%S')
def le_open_left(a, b):
return a is None or (b is not None and a <= b)
def le_open_right(a, b):
return b is None or (a is not None and a <= b)
def eq_open(a, b):
return (a is None and b is None) \
or (a is not None and b is not None and a == b)
def contains(a, b):
return le_open_left(a.start_date, b.start_date) \
and le_open_right(b.end_date, a.end_date)
def find_containing(candidates, node):
for candidate in candidates:
if contains(candidate, node):
return candidate
return None
this_year = time.gmtime()[0]
@ -1615,6 +1641,22 @@ def overlaps(a, b):
)
def check_overlaps(node_type_name, codes, nodes):
errors = []
for ia, a in enumerate(nodes):
for b in nodes[ia+1:]:
if overlaps(a, b):
errors.append(
'%s epochs overlap for %s:\n'
' %s - %s\n %s - %s' % (
node_type_name,
'.'.join(codes),
tts(a.start_date), tts(a.end_date),
tts(b.start_date), tts(b.end_date)))
return errors
class Channel(BaseNode):
'''
Equivalent to SEED blockette 52 and parent element for the related the
@ -2260,15 +2302,8 @@ class FDSNStationXML(Object):
errors = []
for nslc, channels in by_nslc.items():
for ia, a in enumerate(channels):
for b in channels[ia+1:]:
if overlaps(a, b):
errors.append(
'Channel epochs overlap for %s:\n'
' %s - %s\n %s - %s' % (
'.'.join(nslc),
tts(a.start_date), tts(a.end_date),
tts(b.start_date), tts(b.end_date)))
errors.extend(check_overlaps('Channel', nslc, channels))
return errors
def check(self):

185
src/squirrel/base.py

@ -88,6 +88,24 @@ def _is_exact(pat):
return not ('*' in pat or '?' in pat or ']' in pat or '[' in pat)
def prefix_tree(tups):
if not tups:
return []
if len(tups[0]) == 1:
return sorted((tup[0], []) for tup in tups)
d = defaultdict(list)
for tup in tups:
d[tup[0]].append(tup[1:])
sub = []
for k in sorted(d.keys()):
sub.append((k, prefix_tree(d[k])))
return sub
class Batch(object):
'''
Batch of waveforms from window-wise data extraction.
@ -730,7 +748,7 @@ class Squirrel(Selection):
self.add_source(catalog.CatalogSource(*args, **kwargs))
def add_dataset(self, ds, check=True, warn_persistent=True):
def add_dataset(self, ds, check=True):
'''
Read dataset description from file and add its contents.
@ -752,17 +770,6 @@ class Squirrel(Selection):
'''
if isinstance(ds, str):
ds = dataset.read_dataset(ds)
path = ds
else:
path = None
if warn_persistent and ds.persistent and (
not self._persistent or (self._persistent != ds.persistent)):
logger.warning(
'Dataset `persistent` flag ignored. Can not be set on already '
'existing Squirrel instance.%s' % (
' Dataset: %s' % path if path else ''))
ds.setup(self, check=check)
@ -1597,6 +1604,7 @@ class Squirrel(Selection):
try:
return content_cache.get(nut, accessor_id, model)
except KeyError:
raise error.NotAvailable(
'Unable to retrieve content: %s, %s, %s, %s' % nut.key)
@ -1666,16 +1674,6 @@ class Squirrel(Selection):
def get_cache_stats(self, cache_id):
return self._content_caches[cache_id].get_stats()
def _check_duplicates(self, nuts):
d = defaultdict(list)
for nut in nuts:
d[nut.codes].append(nut)
for codes, group in d.items():
if len(group) > 1:
logger.warning(
'Multiple entries matching codes: %s' % str(codes))
@filldocs
def get_stations(
self, obj=None, tmin=None, tmax=None, time=None, codes=None,
@ -1704,13 +1702,13 @@ class Squirrel(Selection):
if model == 'pyrocko':
return self._get_pyrocko_stations(obj, tmin, tmax, time, codes)
elif model in ('squirrel', 'stationxml'):
elif model in ('squirrel', 'stationxml', 'stationxml+'):
args = self._get_selection_args(
STATION, obj, tmin, tmax, time, codes)
nuts = sorted(
self.iter_nuts('station', *args), key=lambda nut: nut.dkey)
self._check_duplicates(nuts)
return [self.get_content(nut, model=model) for nut in nuts]
else:
raise ValueError('Invalid station model: %s' % model)
@ -1736,7 +1734,7 @@ class Squirrel(Selection):
nuts = sorted(
self.iter_nuts('channel', *args), key=lambda nut: nut.dkey)
self._check_duplicates(nuts)
return [self.get_content(nut, model=model) for nut in nuts]
@filldocs
@ -1766,7 +1764,7 @@ class Squirrel(Selection):
nuts = sorted(
self.iter_nuts(
'channel', tmin, tmax, codes), key=lambda nut: nut.dkey)
self._check_duplicates(nuts)
return model.Sensor.from_channels(
self.get_content(nut) for nut in nuts)
@ -1791,7 +1789,7 @@ class Squirrel(Selection):
nuts = sorted(
self.iter_nuts('response', *args), key=lambda nut: nut.dkey)
self._check_duplicates(nuts)
return [self.get_content(nut, model=model) for nut in nuts]
@filldocs
@ -1814,8 +1812,13 @@ class Squirrel(Selection):
See :py:meth:`iter_nuts` for details on time span matching.
'''
if model == 'stationxml':
model_ = 'stationxml+'
else:
model_ = model
responses = self.get_responses(
obj, tmin, tmax, time, codes, model=model)
obj, tmin, tmax, time, codes, model=model_)
if len(responses) == 0:
raise error.NotAvailable(
'No instrument response available (%s).'
@ -1823,11 +1826,15 @@ class Squirrel(Selection):
RESPONSE, obj, tmin, tmax, time, codes))
elif len(responses) > 1:
if model == 'squirrel':
rinfo = ':\n' + '\n'.join(
' ' + resp.summary for resp in responses)
if model_ == 'squirrel':
resps_sq = responses
elif model_ == 'stationxml+':
resps_sq = [resp[0] for resp in responses]
else:
rinfo = '.'
raise ValueError('Invalid response model: %s' % model)
rinfo = ':\n' + '\n'.join(
' ' + resp.summary for resp in resps_sq)
raise error.NotAvailable(
'Multiple instrument responses matching given constraints '
@ -1835,7 +1842,10 @@ class Squirrel(Selection):
self._get_selection_args_str(
RESPONSE, obj, tmin, tmax, time, codes), rinfo))
return responses[0]
if model == 'stationxml':
return responses[0][1]
else:
return responses[0]
@filldocs
def get_events(
@ -1855,7 +1865,7 @@ class Squirrel(Selection):
args = self._get_selection_args(EVENT, obj, tmin, tmax, time, codes)
nuts = sorted(
self.iter_nuts('event', *args), key=lambda nut: nut.dkey)
self._check_duplicates(nuts)
return [self.get_content(nut) for nut in nuts]
def _redeem_promises(self, *args):
@ -2643,6 +2653,113 @@ class Squirrel(Selection):
return coverages
def get_stationxml(
self, obj=None, tmin=None, tmax=None, time=None, codes=None,
level='response'):
'''
Get station/channel/response metadata in StationXML representation.
%(query_args)s
:returns:
:py:class:`~pyrocko.io.stationxml.FDSNStationXML` object.
'''
if level not in ('network', 'station', 'channel', 'response'):
raise ValueError('Invalid level: %s' % level)
tmin, tmax, codes = self._get_selection_args(
CHANNEL, obj, tmin, tmax, time, codes)
filtering = CodesPatternFiltering(codes=codes)
nslcs = list(set(
codes.nslc for codes in
filtering.filter(self.get_codes(kind='channel'))))
from pyrocko.io import stationxml as sx
networks = []
for net, stas in prefix_tree(nslcs):
network = sx.Network(code=net)
networks.append(network)
if level not in ('station', 'channel', 'response'):
continue
for sta, locs in stas:
stations = self.get_stations(
tmin=tmin,
tmax=tmax,
codes=(net, sta, '*'),
model='stationxml')
errors = sx.check_overlaps(
'Station', (net, sta), stations)
if errors:
raise sx.Inconsistencies(
'Inconsistencies found:\n %s'
% '\n '.join(errors))
network.station_list.extend(stations)
if level not in ('channel', 'response'):
continue
for loc, chas in locs:
for cha, _ in chas:
channels = self.get_channels(
tmin=tmin,
tmax=tmax,
codes=(net, sta, loc, cha),
model='stationxml')
errors = sx.check_overlaps(
'Channel', (net, sta, loc, cha), channels)
if errors:
raise sx.Inconsistencies(
'Inconsistencies found:\n %s'
% '\n '.join(errors))
for channel in channels:
station = sx.find_containing(stations, channel)
if station is not None:
station.channel_list.append(channel)
else:
raise sx.Inconsistencies(
'No station or station epoch found for '
'channel: %s' % '.'.join(
(net, sta, loc, cha)))
if level != 'response':
continue
response_sq, response_sx = self.get_response(
codes=(net, sta, loc, cha),
tmin=channel.start_date,
tmax=channel.end_date,
model='stationxml+')
if not (
sx.eq_open(
channel.start_date, response_sq.tmin)
and sx.eq_open(
channel.end_date, response_sq.tmax)):
raise sx.Inconsistencies(
'Response time span does not match '
'channel time span: %s' % '.'.join(
(net, sta, loc, cha)))
channel.response = response_sx
return sx.FDSNStationXML(
source='Generated by Pyrocko Squirrel.',
network_list=networks)
def add_operator(self, op):
self._operators.append(op)

7
src/squirrel/cache.py

@ -125,11 +125,14 @@ class ContentCache(object):
self._accessor_ticks[accessor] = 0
entry[2][accessor] = self._accessor_ticks[accessor]
el = entry[1][element]
if model == 'squirrel':
return entry[1][element].content
return el.content
elif model.endswith('+'):
return el.content, el.raw_content[model[:-1]]
else:
return entry[1][element].raw_content[model]
return el.raw_content[model]
def has(self, nut):
'''

18
src/squirrel/model.py

@ -524,7 +524,7 @@ g_offset_time_unit = 1.0 / g_offset_time_unit_inv
def tsplit(t):
if t is None:
return None, 0.0
return None, 0
t = util.to_time_float(t)
if type(t) is float:
@ -1072,6 +1072,7 @@ class Nut(Object):
self.deltat) = values_nocheck
self.codes = to_codes_simple(self.kind_id, codes_safe_str)
self.deltat = self.deltat or None
self.raw_content = {}
self.content = None
else:
@ -1140,6 +1141,21 @@ class Nut(Object):
self.tmin_seconds, self.tmin_offset,
self.tmax_seconds, self.tmax_offset, self.deltat)
def diff(self, other):
names = [
'file_segment', 'file_element', 'kind_id', 'codes',
'tmin_seconds', 'tmin_offset', 'tmax_seconds', 'tmax_offset',
'deltat']
d = []
for name, a, b in zip(
names, self.equality_values, other.equality_values):
if a != b:
d.append((name, a, b))
return d
@property
def tmin(self):
return tjoin(self.tmin_seconds, self.tmin_offset)

6
src/squirrel/tool/commands/coverage.py

@ -5,6 +5,8 @@
from __future__ import absolute_import, print_function
import time
from pyrocko import util
from pyrocko.plot import terminal
from pyrocko.get_terminal_size import get_terminal_size
@ -36,6 +38,10 @@ def run(parser, args):
tmin_g, tmax_g = squirrel.get_time_span()
sx, _ = get_terminal_size()
now = time.time()
if tmax_g is None or tmax_g > now:
tmax_g = now
kwargs = args.squirrel_query
kinds = kwargs.pop('kind', sq.supported_content_kinds())
codes = kwargs.pop('codes', None)

23
src/squirrel/tool/commands/jackseis.py

@ -112,6 +112,7 @@ class Converter(HasPaths):
out_mseed_steim = IntChoice.T(
optional=True,
choices=[1, 2])
out_meta_path = Path.T(optional=True)
parts = List.T(Defer('Converter.T'))
@ -188,6 +189,12 @@ class Converter(HasPaths):
'Default is STEIM-2, which can compress full range int32. '
'Note: STEIM-2 is limited to 30 bit dynamic range.')
p.add_argument(
'--out-meta-path',
dest='out_meta_path',
metavar='PATH',
help='Set output path for station metadata (StationXML) export.')
@classmethod
def from_arguments(cls, args):
kwargs = args.squirrel_query
@ -200,6 +207,7 @@ class Converter(HasPaths):
out_data_type=args.out_data_type,
out_mseed_record_length=args.out_mseed_record_length,
out_mseed_steim=args.out_mseed_steim,
out_meta_path=args.out_meta_path,
**kwargs)
obj.validate()
@ -264,6 +272,12 @@ class Converter(HasPaths):
else:
return None
def get_effective_out_meta_path(self):
if self.out_meta_path is not None:
return self.expand_path(self.out_meta_path)
else:
return None
def do_rename(self, rules, tr):
rename = {}
for k in ['network', 'station', 'location', 'channel']:
@ -327,6 +341,15 @@ class Converter(HasPaths):
out_format = chain.get('out_format')
out_data_type = chain.get('out_data_type')
out_meta_path = chain.fcall('get_effective_out_meta_path')
if out_meta_path is not None:
sx = sq.get_stationxml(codes=codes, tmin=tmin, tmax=tmax)
util.ensuredirs(out_meta_path)
sx.dump_xml(filename=out_meta_path)
if out_path is None:
return
target_deltat = None
if downsample is not None:
target_deltat = 1.0 / float(downsample)

60
src/squirrel/tool/commands/stationxml.py

@ -0,0 +1,60 @@
# http://pyrocko.org - GPLv3
#
# The Pyrocko Developers, 21st Century
# ---|P------/S----------~Lg----------
from __future__ import absolute_import, print_function
from pyrocko.squirrel.tool.common import ldq, dq
headline = 'Export station metadata as StationXML.'
def make_subparser(subparsers):
return subparsers.add_parser(
'stationxml',
help=headline,
description=headline)
def setup(parser):
parser.add_squirrel_selection_arguments()
parser.add_squirrel_query_arguments()
level_choices = ('network', 'station', 'channel', 'response')
level_default = 'response'
parser.add_argument(
'--level',
dest='level',
choices=level_choices,
default=level_default,
help='Set level of detail to be returned. Choices: %s. '
'Default: %s.' % (
ldq(level_choices),
dq(level_default)))
out_format_choices = ('xml', 'yaml')
out_format_default = 'xml'
parser.add_argument(
'--out-format',
dest='out_format',
choices=out_format_choices,
default=out_format_default,
help='Set format of output. Choices: %s. '
'Default: %s.' % (
ldq(out_format_choices),
dq(out_format_default)))
def run(parser, args):
squirrel = args.make_squirrel()
sx = squirrel.get_stationxml(**args.squirrel_query, level=args.level)
if args.out_format == 'xml':
print(sx.dump_xml())
elif args.out_format == 'yaml':
print(sx.dump())
else:
assert False
Loading…
Cancel
Save