Dataset management - Squirrel

At the base of Pyrocko’s waveform dataset handling is the pyrocko.squirrel framework, centered around the Squirrel class. It indexes and caches the meta-data of local and remote waveform datasets, however large or however many files they are split into, and provides on-demand loading of waveform data so that only data relevant to the current process/view has to be read from disk (or downloaded) into the limited computer memory. It can efficiently handle selections with millions of files. See Squirrel - Data access infrastructure for a conceptual overview.

Reorganizing a dataset into hour-files

In each iteration we get all data for the current time window as a Batch, yielded by pyrocko.squirrel.base.Squirrel.chopper_waveforms(). The batch ‘knows’ the time window to which it belongs, stored in the attributes batch.tmin and batch.tmax (note: a trace’s own onset, trace.tmin, does not have to be identical to batch.tmin). The directory parts in the output path will be created as necessary. When applying this procedure to a dataset consisting of arbitrarily separated files, it will automatically connect adjacent traces as needed!

from pyrocko import io, util
from pyrocko.squirrel import Squirrel
from pyrocko.example import get_example_data

''' Chop a dataset into one-hour windows and reorganize it into hour-files.

For simple cases like this, the same result can also be achieved without
writing any code, using the command line tool `squirrel jackseis`:

    squirrel jackseis --add test.mseed --tinc 3600 \\
        --out-path 'test_hourfiles/hourfile-%(wmin)s.mseed'
'''

# Download test file
get_example_data('test.mseed')

# All data access happens through a Squirrel instance. Could give
# directories or thousands of filenames here, or even add online data
# sources with sq.add_fdsn(...).
sq = Squirrel()
sq.add(['test.mseed'])

# Get timestamp for full hour before first data sample in all selected
# traces.
tmin, _ = sq.get_time_span()
tmin = util.hour_start(tmin)

# Iterate over the data, with a window length of one hour. This will
# automatically connect adjacent traces from separate files as needed.
for batch in sq.chopper_waveforms(tmin=tmin, tinc=3600.):
    if batch.traces:  # the list could be empty due to gaps
        timestring = util.time_to_str(batch.tmin, format='%Y-%m-%d_%H')
        filepath = 'test_hourfiles/hourfile-%s.mseed' % timestring
        io.save(batch.traces, filepath)

Download squirrel_hour_files.py

Downsampling a whole dataset

Example for downsampling all trace files in the input folder to a common sampling rate with pyrocko.trace.Trace.downsample_to(). Padding (tpad) is added around each window and trimmed off again after downsampling, to absorb filter edge effects at window boundaries.

from pyrocko import io, trace, util
from pyrocko.squirrel import Squirrel
from pyrocko.example import get_example_data

'''
Downsample a whole dataset to a common sampling rate.

For simple cases like this, the same result can also be achieved without
writing any code, using the command line tool `squirrel jackseis`:

    squirrel jackseis --add test.mseed --downsample 0.5 --tinc 3600 \\
        --out-path 'downsampled/%(station)s_%(channel)s_%(wmin)s.mseed'
'''

# Download test file
get_example_data('test.mseed')

sq = Squirrel()
sq.add(['test.mseed'])

tinc = 3600.
target_deltat = 2.0  # test.mseed is sampled at 1 Hz (deltat=1s)
tpad = 50 * target_deltat

# Iterate over the data, with a window length of one hour and padding on
# either side to absorb downsampling filter edge effects. Start time windows
# at full hours.
for batch in sq.chopper_waveforms(tinc=tinc, tpad=tpad, snap_window=True):
    traces = []
    for tr in batch.traces:
        tr.downsample_to(target_deltat, snap=True, demean=False)

        try:
            # remove padding
            tr.chop(batch.tmin, batch.tmax)
            traces.append(tr)
        except trace.NoData:
            # can happen for the trailing, incomplete window
            pass

    if traces:  # the list could be empty due to gaps
        timestring = util.time_to_str(batch.tmin, format='%Y-%m-%d_%H')
        filepath = 'downsampled/%(station)s_%(channel)s_%(mytimestring)s.mseed'
        io.save(traces, filepath, additional={'mytimestring': timestring})


# now look at the result with
#   > squirrel snuffler --add downsampled/

Download squirrel_downsample.py

Converting a dataset from Mini-SEED to SAC format

Plain format conversion, possibly combined with time-windowing, downsampling, or restitution, does not require writing any Python code - the command line tool squirrel jackseis handles this directly. See pyrocko.io for the list of formats supported for output.

squirrel jackseis --add mseed/ --tinc 3600 --out-format sac \
    --out-path 'sac/%(station)s_%(channel)s_%(wmin)s.sac'

If station meta-data (for example to fill in station coordinates and channel orientation into the SAC headers) needs to be attached while converting with pyrocko.io.save() directly, pass a dict of Station objects, keyed by (network, station, location), as its stations argument.