mth5.timeseries package

Subpackages

Submodules

mth5.timeseries.channel_ts module

Channel time series module for MT data.

This module provides the ChannelTS class for handling magnetotelluric (MT) time series data with comprehensive metadata management, calibration, and signal processing capabilities.

Notes

  • Time series are stored in xarray.DataArray for efficient operations.

  • Metadata follows the mt_metadata standard with Survey/Station/Run/Channel hierarchy.

  • Supports instrument response removal, resampling, merging, and Obspy integration.

class mth5.timeseries.channel_ts.ChannelTS(channel_type: str = 'auxiliary', data: ndarray | DataFrame | Series | DataArray | list | tuple | None = None, channel_metadata: Electric | Magnetic | Auxiliary | dict | None = None, station_metadata: Station | dict | None = None, run_metadata: Run | dict | None = None, survey_metadata: Survey | dict | None = None, **kwargs: Any)[source]

Bases: object

Time series container for a single MT channel with full metadata.

Stores equally-spaced time series data in an xarray.DataArray with a time coordinate index. Integrates comprehensive metadata from Survey/Station/Run/Channel hierarchy and supports calibration, resampling, merging, and format conversions.

Parameters:
  • channel_type ({'electric', 'magnetic', 'auxiliary'}, default 'auxiliary') – Type of the channel.

  • data (array-like, optional) – Time series data (numpy array, pandas DataFrame/Series, xarray.DataArray).

  • channel_metadata (mt_metadata.timeseries.Electric | Magnetic | Auxiliary | dict, optional) – Channel-specific metadata.

  • station_metadata (mt_metadata.timeseries.Station | dict, optional) – Station metadata.

  • run_metadata (mt_metadata.timeseries.Run | dict, optional) – Run metadata.

  • survey_metadata (mt_metadata.timeseries.Survey | dict, optional) – Survey metadata.

  • **kwargs – Additional attributes to set on the object.

ts[source]

The time series data array.

Type:

numpy.ndarray

sample_rate[source]

Sample rate in samples per second.

Type:

float

start[source]

Start time (UTC).

Type:

MTime

end[source]

End time (UTC), derived from start + duration.

Type:

MTime

n_samples[source]

Number of samples.

Type:

int

component[source]

Component name (e.g., ‘ex’, ‘hy’, ‘temperature’).

Type:

str

channel_response[source]

Full instrument response filter chain.

Type:

ChannelResponse

Notes

  • End time is a derived property and cannot be set directly.

  • Leverages xarray for efficient interpolation, resampling, and groupby operations.

  • Metadata follows mt_metadata standards with automatic time period updates.

Examples

Create an auxiliary channel with synthetic data:

>>> from mth5.timeseries import ChannelTS
>>> import numpy as np
>>> ts_obj = ChannelTS('auxiliary')
>>> ts_obj.sample_rate = 8
>>> ts_obj.start = '2020-01-01T12:00:00+00:00'
>>> ts_obj.ts = np.random.randn(4096)
>>> ts_obj.station_metadata.id = 'MT001'
>>> ts_obj.run_metadata.id = 'MT001a'
>>> ts_obj.component = 'temperature'
>>> print(ts_obj)

Calibrate and remove instrument response:

>>> calibrated = ts_obj.remove_instrument_response()
>>> calibrated.channel_metadata.units
property channel_metadata: Electric | Magnetic | Auxiliary[source]

Channel metadata.

Returns:

Channel metadata from the first channel in the run.

Return type:

mt_metadata.timeseries.Electric | Magnetic | Auxiliary

property channel_response[source]

Full channel response filter

Returns:

full channel response filter

Return type:

mt_metadata.timeseries.filters.ChannelResponse

property channel_type: str[source]

Channel type.

Returns:

Channel type: ‘Electric’, ‘Magnetic’, or ‘Auxiliary’.

Return type:

str

property component[source]
compute_sample_rate()[source]

Two cases, high_frequency (HF) data and not HF data.

# Original comment about the HF case: Taking the median(diff(timestamps)) is more accurate for high sample rates, the way pandas.date_range rounds nanoseconds is not consistent between samples, therefore taking the median provides better results if the time series is long this can be inefficient so test first

copy(data: bool = True) ChannelTS[source]

Create a copy of the ChannelTS object.

Parameters:

data (bool, default True) – Include data in the copy (True) or only metadata (False).

Returns:

Copy of the channel.

Return type:

ChannelTS

Examples

Copy metadata structure without data:

>>> ch_copy = ts_obj.copy(data=False)
decimate(new_sample_rate, inplace=False, max_decimation=8)[source]

decimate the data by using scipy.signal.decimate

Parameters:

dec_factor (int) – decimation factor

  • refills ts.data with decimated data and replaces sample_rate

property end[source]

MTime object

from_obspy_trace(obspy_trace)[source]

Fill data from an obspy.core.Trace

Parameters:

obspy_trace (obspy.core.trace) – Obspy trace object

get_calibrated_units()[source]

Follows the FDSN standard which has the filter stages starting with physical units to digital counts.

The channel_response is expected to have a list of filter “stages” of which the first stage has input units corresponding to the the physical quantity that the instrument measures, and the last is normally counts.

channel_response can be viewed as the chaining together of all of these filters.

Thus it is normal for channel_response.units_out will be in the same units as the archived raw time series, and for the units after the response is corrected for will be the units_in of

The units of the channel metadata are compared to the input and output units of the channel_response.

Returns:

tuple, calibration_operation, either “mulitply” or divide”, and a string for calibrated units

Return type:

tuple (of two strings_

get_calibration_operation()[source]
get_sample_rate_supplied_at_init(channel_metadata: Electric | Magnetic | Auxiliary | dict | None) float | None[source]

Extract sample_rate from channel_metadata if available.

Parameters:

channel_metadata (mt_metadata.timeseries.Electric | Magnetic | Auxiliary | dict | None) – Metadata that may contain a sample_rate field.

Returns:

Sample rate if found, otherwise None.

Return type:

float | None

Notes

Supports nested dict structures like {"electric": {"sample_rate": 8.0}}.

get_slice(start, end=None, n_samples=None)[source]

Get a slice from the time series given a start and end time.

Looks for >= start & <= end

Uses loc to be exact with milliseconds

Parameters:
  • start (string, MTime) – start time of the slice

  • end (string, MTime) – end time of the slice

  • n_samples (integer) – number of sample to get after start time

Returns:

slice of the channel requested

Return type:

ChannelTS

has_data()[source]

check to see if there is an index in the time series

is_high_frequency(threshold_dt=0.0001)[source]

Quasi hard-coded condition to check if data are logged at more than 10kHz can be parameterized in future

merge(other, gap_method='slinear', new_sample_rate=None, resample_method='poly')[source]

merg two channels or list of channels together in the following steps

  1. xr.combine_by_coords([original, other])

  2. compute monotonic time index

  3. reindex(new_time_index, method=gap_method)

If you want a different method or more control use merge

Parameters:

other (mth5.timeseries.ChannelTS) – Another channel

Raises:
  • TypeError – If input is not a ChannelTS

  • ValueError – if the components are different

Returns:

Combined channel with monotonic time index and same metadata

Return type:

mth5.timeseries.ChannelTS

property n_samples[source]

number of samples

plot()[source]

Simple plot of the data

Returns:

figure object

Return type:

matplotlib.figure

plot_spectra(spectra_type='welch', window_length=4096, **kwargs)[source]
Parameters:
  • spectra_type (string, optional) – spectra type, defaults to “welch”

  • window_length (int, optional) – window length of the welch method should be a power of 2, defaults to 2 ** 12

  • **kwargs

    DESCRIPTION

remove_instrument_response(include_decimation=False, include_delay=False, **kwargs)[source]

Remove instrument response from the given channel response filter

The order of operations is important (if applied):

  1. detrend

  2. zero mean

  3. zero pad

  4. time window

  5. frequency window

  6. remove response

  7. undo time window

  8. bandpass

Parameters:
  • include_decimation (bool, optional) – Include decimation in response, defaults to True

  • include_delay (bool, optional) – include delay in complex response, defaults to False

kwargs

Parameters:
  • plot (boolean, default True) – to plot the calibration process [ False | True ]

  • detrend (boolean, default True) – Remove linar trend of the time series

  • zero_mean (boolean, default True) – Remove the mean of the time series

  • zero_pad (boolean, default True) – pad the time series to the next power of 2 for efficiency

  • t_window (string, default None) – Time domain windown name see scipy.signal.windows for options

  • t_window_params – Time domain window parameters, parameters can be

found in scipy.signal.windows :type t_window_params: dictionary :param f_window: Frequency domain windown name see scipy.signal.windows for options :type f_window: string, defualt None :param f_window_params: Frequency window parameters, parameters can be found in scipy.signal.windows :type f_window_params: dictionary :param bandpass: bandpass freequency and order {“low”:, “high”:, “order”:,} :type bandpass: dictionary

resample_poly(new_sample_rate, pad_type='mean', inplace=False)[source]

Use scipy.signal.resample_poly to resample data while using an FIR filter to remove aliasing.

Parameters:
  • new_sample_rate (TYPE) – DESCRIPTION

  • pad_type (TYPE, optional) – DESCRIPTION, defaults to “mean”

Returns:

DESCRIPTION

Return type:

TYPE

property run_metadata: Run[source]

Run metadata.

Returns:

Run metadata from the first run in the station.

Return type:

mt_metadata.timeseries.Run

property sample_interval[source]

Sample interval = 1 / sample_rate

Returns:

sample interval as time distance between time samples

Return type:

float

property sample_rate[source]

sample rate in samples/second

property start[source]

MTime object

property station_metadata: Station[source]

Station metadata.

Returns:

Station metadata from the first station in the survey.

Return type:

mt_metadata.timeseries.Station

property survey_metadata: Survey[source]

Survey metadata.

Returns:

Survey metadata with updated keys.

Return type:

mt_metadata.timeseries.Survey

property time_index: ndarray[source]

Time index as a numpy array.

Returns:

Array of datetime64[ns] timestamps.

Return type:

numpy.ndarray

to_obspy_trace(network_code=None, encoding=None)[source]

Convert the time series to an obspy.core.trace.Trace object. This will be helpful for converting between data pulled from IRIS and data going into IRIS.

Parameters:

network_code (string) – two letter code provided by FDSN DMC

Returns:

DESCRIPTION

Return type:

TYPE

to_xarray()[source]

Returns a xarray.DataArray object of the channel timeseries this way metadata from the metadata class is updated upon return.

Returns:

Returns a xarray.DataArray object of the channel timeseries

this way metadata from the metadata class is updated upon return. :rtype: xarray.DataArray

>>> import numpy as np
>>> from mth5.timeseries import ChannelTS
>>> ex = ChannelTS("electric")
>>> ex.start = "2020-01-01T12:00:00"
>>> ex.sample_rate = 16
>>> ex.ts = np.random.rand(4096)
property ts: ndarray[source]

Time series data as a numpy array.

Returns:

The time series data.

Return type:

numpy.ndarray

welch_spectra(window_length=4096, **kwargs)[source]

get welch spectra

Parameters:
  • window_length (TYPE) – DESCRIPTION

  • **kwargs

    DESCRIPTION

Returns:

DESCRIPTION

Return type:

TYPE

mth5.timeseries.run_ts module

lists and arrays that goes on, seems easiest to convert all lists to strings and then convert them back if read in.

copyright:

Jared Peacock (jpeacock@usgs.gov)

license:

MIT

class mth5.timeseries.run_ts.RunTS(array_list: list[ChannelTS] | list[DataArray] | Dataset | None = None, run_metadata: Run | dict | None = None, station_metadata: Station | dict | None = None, survey_metadata: Survey | dict | None = None)[source]

Bases: object

Container for MT time series data from a single run.

Holds all run time series in one aligned xarray Dataset with channels as data variables and time as the coordinate. Manages metadata for survey, station, and run levels.

Parameters:
  • array_list (list[ChannelTS] | list[xr.DataArray] | xr.Dataset | None, optional) – List of ChannelTS objects, xarray DataArrays, or an xarray Dataset containing the time series data. All channels will be aligned to a common time index.

  • run_metadata (timeseries.Run | dict | None, optional) – Metadata for the run. Can be a Run object or dictionary.

  • station_metadata (timeseries.Station | dict | None, optional) – Metadata for the station. Can be a Station object or dictionary.

  • survey_metadata (timeseries.Survey | dict | None, optional) – Metadata for the survey. Can be a Survey object or dictionary.

dataset[source]

xarray Dataset containing all channel data with time coordinate

Type:

xr.Dataset

survey_metadata[source]

Survey-level metadata

Type:

timeseries.Survey

station_metadata[source]

Station-level metadata

Type:

timeseries.Station

run_metadata[source]

Run-level metadata

Type:

timeseries.Run

filters[source]

Dictionary of channel response filters keyed by filter name

Type:

dict[str, Filter]

sample_rate[source]

Sample rate in samples per second

Type:

float

channels[source]

List of channel names in the dataset

Type:

list[str]

Examples

Create an empty RunTS:

>>> from mth5.timeseries import RunTS
>>> run = RunTS()

Create RunTS from ChannelTS objects:

>>> from mth5.timeseries import ChannelTS, RunTS
>>> ex = ChannelTS('electric', data=ex_data,
...                channel_metadata={'component': 'ex'})
>>> ey = ChannelTS('electric', data=ey_data,
...                channel_metadata={'component': 'ey'})
>>> run = RunTS(array_list=[ex, ey])
>>> print(run.channels)
['ex', 'ey']

Access individual channels:

>>> ex_channel = run.ex  # Returns ChannelTS object
>>> print(ex_channel.sample_rate)
256.0

See also

ChannelTS

Individual channel time series container

Notes

When multiple channels are provided with different start/end times, they will be automatically aligned using the earliest start and latest end times, with NaN values filling gaps.

add_channel(channel: DataArray | ChannelTS) None[source]

Add a channel to the dataset.

The channel must have the same sample rate and time coordinates that are compatible with the existing dataset. If start times don’t match, NaN values will be placed where timing doesn’t align.

Parameters:

channel (xr.DataArray | ChannelTS) – A channel as an xarray DataArray or ChannelTS object to add.

Raises:

ValueError – If the channel has a different sample rate than the run, or if the input is not a DataArray or ChannelTS.

Examples

Add a ChannelTS:

>>> hz = ChannelTS('magnetic', data=hz_data,
...                channel_metadata={'component': 'hz'})
>>> run.add_channel(hz)
>>> print(run.channels)
['ex', 'ey', 'hx', 'hy', 'hz']

Add an xarray DataArray:

>>> import xarray as xr
>>> data_array = xr.DataArray(data, coords={'time': times})
>>> run.add_channel(data_array)
calibrate(**kwargs) RunTS[source]

Remove instrument response from all channels.

Applies the channel response filters to calibrate each channel, creating a new run with calibrated data.

Parameters:

**kwargs – Additional keyword arguments passed to each channel’s remove_instrument_response method.

Returns:

New RunTS object with calibrated channels.

Return type:

RunTS

Examples

>>> calibrated_run = run.calibrate()
>>> # Calibration typically converts from counts to physical units

See also

ChannelTS.remove_instrument_response

Calibrate single channel

property channels: list[str][source]

List of channel names in the dataset.

Returns:

List of channel component names (e.g., [‘ex’, ‘ey’, ‘hx’]).

Return type:

list[str]

Examples

>>> print(run.channels)
['ex', 'ey', 'hx', 'hy', 'hz']
copy(data: bool = True) RunTS[source]

Create a copy of the RunTS object.

Parameters:

data (bool, optional) – If True, copy the data along with timeseries. If False, only copy the metadata (default is True).

Returns:

A copy of the RunTS object.

Return type:

RunTS

Examples

Create a copy with data:

>>> run_copy = run.copy()

Create a metadata-only copy:

>>> run_meta = run.copy(data=False)
>>> print(run_meta.has_data())
False
property dataset: Dataset[source]

The xarray Dataset containing all channel data.

Returns:

Dataset with channels as data variables and time as coordinate.

Return type:

xr.Dataset

Examples

>>> print(run.dataset)
<xarray.Dataset>
Dimensions:  (time: 4096)
Coordinates:
  * time     (time) datetime64[ns] ...
Data variables:
    ex       (time) float64 ...
    ey       (time) float64 ...
decimate(new_sample_rate: float, inplace: bool = False, max_decimation: int = 8) RunTS | None[source]

Decimate data to a new sample rate using multi-stage decimation.

Applies FIR filtering and downsampling in multiple stages to achieve the target sample rate while preventing aliasing.

Parameters:
  • new_sample_rate (float) – Target sample rate in samples per second.

  • inplace (bool, optional) – If True, modify the current run. If False, return a new run (default is False).

  • max_decimation (int, optional) – Maximum decimation factor for each stage (default is 8).

Returns:

If inplace is False, returns new decimated RunTS. Otherwise None.

Return type:

RunTS | None

Examples

Decimate from 256 Hz to 1 Hz:

>>> decimated_run = run.decimate(1.0)
>>> print(decimated_run.sample_rate)
1.0

Decimate in place:

>>> run.decimate(16.0, inplace=True)
>>> print(run.sample_rate)
16.0

Notes

NaN values are filled with 0 before decimation to prevent NaN propagation. Multi-stage decimation is used to maintain signal quality and prevent aliasing.

See also

resample_poly

Alternative resampling using polyphase filtering

resample

Simple resampling without anti-aliasing

property end: MTime[source]

End time of the run in UTC.

Returns:

End time from the dataset if data exists, otherwise from run_metadata.

Return type:

MTime

Examples

>>> print(run.end)
2020-01-01T01:00:00+00:00
property filters: dict[str, ChannelResponse][source]

Dictionary of channel response filters.

Returns:

Dictionary keyed by filter name containing ChannelResponse objects.

Return type:

dict[str, ChannelResponse]

Examples

>>> print(run.filters.keys())
dict_keys(['v_to_counts', 'dipole_100m'])
from_obspy_stream(obspy_stream: Stream, run_metadata: Run | None = None) None[source]

Get a run from an obspy.core.stream which is a list of obspy.core.Trace objects.

Parameters:

obspy_stream (obspy.core.Stream) – Obspy Stream object

Development Notes:
  • There is a baked in assumption here that the channel nomenclature in obspy is e1,e2,h1,h2,h3 and we want to convert to mth5 conventions ex,ey,hx,hy,hz. This should be made more flexible in the future.

  • A bug was found that was creating channels e1, ex, ey in the same run when reading from obspy – this is fixed here by renaming the components and a workaround to reset the station’s channels_recorded list.

get_slice(start: str | MTime, end: str | MTime | None = None, n_samples: int | None = None) RunTS[source]

Extract a time slice from the run.

Gets a chunk of data from the run, finding the closest points to the given parameters. Uses pandas slice_indexer for robust slicing.

Parameters:
  • start (str | MTime) – Start time of the slice (ISO format string or MTime object).

  • end (str | MTime | None, optional) – End time of the slice. Required if n_samples not provided.

  • n_samples (int | None, optional) – Number of samples to get. Required if end not provided.

Returns:

New RunTS object containing the requested slice with copies of metadata and filters.

Return type:

RunTS

Raises:

ValueError – If neither end nor n_samples is provided.

Examples

Get slice by start and end times:

>>> slice1 = run.get_slice('2020-01-01T00:00:00',
...                         '2020-01-01T00:01:00')
>>> print(slice1.start, slice1.end)

Get slice by start time and number of samples:

>>> slice2 = run.get_slice('2020-01-01T00:00:00', n_samples=1024)
>>> print(len(slice2.dataset.time))
1024

Notes

Uses pandas slice_indexer which handles near-matches better than xarray’s native slicing. The actual slice may be slightly adjusted to match available data points.

has_data() bool[source]

Check if the RunTS contains any data.

Returns:

True if channels with data exist, False otherwise.

Return type:

bool

Examples

>>> run = RunTS()
>>> print(run.has_data())
False
>>> run.add_channel(ex_channel)
>>> print(run.has_data())
True
merge(other: RunTS | list[RunTS], gap_method: str = 'slinear', new_sample_rate: float | None = None, resample_method: str = 'poly') RunTS[source]

Merge multiple runs into a single run.

Combines this run with one or more other runs, optionally resampling to a common sample rate and filling gaps with interpolation.

Parameters:
  • other (RunTS | list[RunTS]) – Another RunTS object or list of RunTS objects to merge.

  • gap_method (str, optional) – Interpolation method for filling gaps: ‘linear’, ‘nearest’, ‘zero’, ‘slinear’, ‘quadratic’, ‘cubic’ (default is ‘slinear’).

  • new_sample_rate (float | None, optional) – If provided, all runs will be resampled to this rate before merging. If None, uses the sample rate of the first run.

  • resample_method (str, optional) – Resampling method if new_sample_rate is provided: ‘decimate’ or ‘poly’ (default is ‘poly’).

Returns:

New merged RunTS object with monotonic time index.

Return type:

RunTS

Raises:

TypeError – If other is not a RunTS or list of RunTS objects.

Examples

Merge two runs:

>>> run1 = RunTS(array_list=[ex1, ey1])
>>> run2 = RunTS(array_list=[ex2, ey2])
>>> merged = run1.merge(run2)

Merge multiple runs with resampling:

>>> runs = [run1, run2, run3]
>>> merged = run1.merge(runs, new_sample_rate=1.0,
...                     resample_method='poly')

Notes

The merge process:

  1. Optionally resample all runs to common sample rate

  2. Combine datasets using xr.combine_by_coords

  3. Create monotonic time index spanning full range

  4. Interpolate to new index filling gaps

  5. Merge all filter dictionaries

Metadata is taken from the first run (self).

See also

__add__

Simple merging with + operator

plot(color_map: dict[str, tuple[float, float, float]] | None = None, channel_order: list[str] | None = None) Figure[source]

Plot all channels as time series.

Creates a multi-panel figure with each channel in its own subplot, sharing a common time axis.

Parameters:
  • color_map (dict[str, tuple[float, float, float]] | None, optional) –

    Dictionary mapping channel names to RGB color tuples (values 0-1). Default colors:

    • ex: (1, 0.2, 0.2) - red

    • ey: (1, 0.5, 0) - orange

    • hx: (0, 0.5, 1) - blue

    • hy: (0.5, 0.2, 1) - purple

    • hz: (0.2, 1, 1) - cyan

  • channel_order (list[str] | None, optional) – Order of channels from top to bottom. If None, uses order from self.channels.

Returns:

Figure object containing the plot.

Return type:

matplotlib.figure.Figure

Examples

Plot with default settings:

>>> fig = run.plot()
>>> fig.savefig('timeseries.png')

Plot with custom colors and order:

>>> colors = {'ex': (1, 0, 0), 'ey': (0, 1, 0)}
>>> fig = run.plot(color_map=colors, channel_order=['ey', 'ex'])

Warning

May be slow for large datasets (millions of points). Consider using get_slice() first to plot a subset.

plot_spectra(spectra_type: str = 'welch', color_map: dict[str, tuple[float, float, float]] | None = None, **kwargs) Figure[source]

Plot power spectral density for all channels.

Computes and plots the power spectrum of each channel on a single log-log plot with period on x-axis.

Parameters:
  • spectra_type (str, optional) – Type of spectral estimate to compute. Currently only ‘welch’ is supported (default is ‘welch’).

  • color_map (dict[str, tuple[float, float, float]] | None, optional) – Dictionary mapping channel names to RGB color tuples (values 0-1). Uses same defaults as plot().

  • **kwargs – Additional keyword arguments passed to the spectra computation method (e.g., nperseg, window for Welch’s method).

Returns:

Figure object containing the spectra plot.

Return type:

matplotlib.figure.Figure

Examples

Plot spectra with default settings:

>>> fig = run.plot_spectra()

Plot with custom Welch parameters:

>>> fig = run.plot_spectra(nperseg=1024, window='hann')

Notes

The plot shows:

  • Period (seconds) on bottom x-axis

  • Frequency (Hz) on top x-axis

  • Power (dB) on y-axis

See also

ChannelTS.welch_spectra

Compute Welch power spectrum

resample(new_sample_rate: float, inplace: bool = False) RunTS | None[source]

Resample data to a new sample rate using nearest-neighbor method.

Simple resampling without anti-aliasing filtering. Use decimate or resample_poly for better quality when downsampling.

Parameters:
  • new_sample_rate (float) – Target sample rate in samples per second.

  • inplace (bool, optional) – If True, modify current run. If False, return new run (default is False).

Returns:

If inplace is False, returns new resampled RunTS. Otherwise None.

Return type:

RunTS | None

Examples

>>> resampled_run = run.resample(128.0)
>>> print(resampled_run.sample_rate)
128.0

Warning

This method does not apply anti-aliasing filtering. For downsampling, consider using decimate() or resample_poly() instead.

See also

decimate

Proper downsampling with anti-aliasing

resample_poly

High-quality resampling with polyphase filtering

resample_poly(new_sample_rate: float, pad_type: str = 'mean', inplace: bool = False) RunTS | None[source]

Resample data using polyphase filtering.

Uses scipy.signal.resample_poly to resample while applying an FIR filter to remove aliasing. Generally more accurate than simple resampling but slower than decimation.

Parameters:
  • new_sample_rate (float) – Target sample rate in samples per second.

  • pad_type (str, optional) – Padding method for edge effects: ‘mean’, ‘median’, ‘zero’ (default is ‘mean’).

  • inplace (bool, optional) – If True, modify current run. If False, return new run (default is False).

Returns:

If inplace is False, returns new resampled RunTS. Otherwise None.

Return type:

RunTS | None

Examples

Resample from 256 Hz to 100 Hz:

>>> resampled_run = run.resample_poly(100.0)
>>> print(resampled_run.sample_rate)
100.0

Notes

NaN values are filled with 0 before resampling. The polyphase method is particularly good for arbitrary sample rate ratios.

See also

decimate

Multi-stage decimation for downsampling

resample

Simple nearest-neighbor resampling

property run_metadata: Run[source]

Run timeseries.

Returns:

Run-level metadata object (first run in first station).

Return type:

timeseries.Run

property sample_interval: float[source]

Sample interval in seconds (inverse of sample_rate).

Returns:

Sample interval = 1 / sample_rate, or 0.0 if sample_rate is 0.

Return type:

float

Examples

>>> print(run.sample_interval)
0.00390625  # for 256 Hz
property sample_rate: float[source]

Sample rate in samples per second.

Returns:

Sample rate estimated from time differences if data exists, otherwise from timeseries.

Return type:

float

Examples

>>> print(run.sample_rate)
256.0
set_dataset(array_list: list[ChannelTS] | list[DataArray] | Dataset, align_type: str = 'outer') None[source]

Set the dataset from a list of channels or existing dataset.

Creates an xarray Dataset from the input channels or dataset, validates metadata consistency, and updates dataset attributes with run metadata.

Parameters:
  • array_list (list[ChannelTS] | list[xr.DataArray] | xr.Dataset) – Input data as a list of ChannelTS objects, list of xarray DataArrays, or an existing xarray Dataset.

  • align_type (str, optional) –

    Method for aligning channels with different time indexes:

    • ’outer’ - use union of all time indexes (default)

    • ’inner’ - use intersection of time indexes

    • ’left’ - use time index from first channel

    • ’right’ - use time index from last channel

    • ’exact’ - raise ValueError if indexes don’t match exactly

    • ’override’ - rewrite indexes to match first channel (requires same size)

Notes

This method performs the following operations:

  1. Validates the input array_list

  2. Converts ChannelTS objects to xarray format

  3. Combines channels into a single Dataset

  4. Validates metadata consistency

  5. Updates dataset attributes with run metadata

When providing ChannelTS objects, their metadata and filters are automatically extracted and merged into the run’s metadata structure.

Examples

Set dataset from ChannelTS objects:

>>> ex = ChannelTS('electric', data=ex_data,
...                channel_metadata={'component': 'ex'})
>>> ey = ChannelTS('electric', data=ey_data,
...                channel_metadata={'component': 'ey'})
>>> run.set_dataset([ex, ey])

Set dataset with custom alignment:

>>> run.set_dataset([ex, ey, hx], align_type='inner')

Set dataset from existing xarray Dataset:

>>> import xarray as xr
>>> ds = xr.Dataset({'ex': ex_da, 'ey': ey_da})
>>> run.set_dataset(ds)

See also

dataset

Property for setting dataset with default alignment

add_channel

Add a single channel to existing dataset

_validate_array_list

Validation and conversion of array list

property start: MTime[source]

Start time of the run in UTC.

Returns:

Start time from the dataset if data exists, otherwise from run_metadata.

Return type:

MTime

Examples

>>> print(run.start)
2020-01-01T00:00:00+00:00
property station_metadata: Station[source]

Station timeseries.

Returns:

Station-level metadata object (first station in survey).

Return type:

timeseries.Station

property summarize_metadata: dict[str, any][source]

Get a summary of all channel timeseries.

Flattens the metadata from all channels into a single dictionary with keys in the format ‘channel.attribute’.

Returns:

Dictionary with flattened metadata from all channels.

Return type:

dict[str, any]

Examples

>>> meta_summary = run.summarize_metadata
>>> print(meta_summary.keys())
dict_keys(['ex.time_period.start', 'ex.sample_rate', ...])
property survey_metadata: Survey[source]

Survey timeseries.

Returns:

Survey-level metadata object.

Return type:

timeseries.Survey

to_obspy_stream(network_code: str | None = None, encoding: str | None = None) Stream[source]

Convert time series to an ObsPy Stream object.

Creates an ObsPy Stream containing a Trace for each channel in the run.

Parameters:
  • network_code (str | None, optional) – Two-letter network code provided by FDSN DMC. If None, uses station timeseries.

  • encoding (str | None, optional) – Data encoding format (e.g., ‘STEIM2’, ‘FLOAT64’). If None, uses default encoding.

Returns:

Stream object containing Trace objects for all channels.

Return type:

obspy.core.Stream

Examples

>>> stream = run.to_obspy_stream(network_code='MT')
>>> print(stream)
3 Trace(s) in Stream:
MT.MT001..EX | 2020-01-01T00:00:00 - ... | 256.0 Hz, 4096 samples

See also

from_obspy_stream

Create RunTS from ObsPy Stream

ChannelTS.to_obspy_trace

Convert single channel

validate_metadata() None[source]

Check to make sure that the metadata matches what is in the data set.

updates metadata from the data.

Check the start and end times, channels recorded

wrangle_leap_seconds_from_obspy(array_list: list[ChannelTS]) list[ChannelTS][source]

Handle potential leap second issues from ObsPy streams.

Removes runs with only one sample that are numerically identical to adjacent samples, which may be artifacts of leap second handling.

Parameters:

array_list (list[ChannelTS]) – List of ChannelTS objects from ObsPy conversion.

Returns:

Filtered list with single-sample runs removed.

Return type:

list[ChannelTS]

Notes

This is experimental handling for issue #169. The exact behavior of ObsPy’s leap second handling is not fully documented.

mth5.timeseries.scipy_filters module

Scipy filter wrappers for xarray.

This module provides xarray-compatible wrappers for scipy.signal filtering functions, enabling efficient filtering operations on labeled N-dimensional arrays with automatic dimension handling.

Notes

  • Adapted from xr-scipy: https://github.com/fujiisoup/xr-scipy

  • Filters can be applied along any dimension with automatic sampling rate detection.

  • Supports both IIR and FIR filter types with forward/backward filtering.

Examples

Apply a bandpass filter to an xarray DataArray:

>>> import xarray as xr
>>> import numpy as np
>>> data = xr.DataArray(np.random.randn(1000), dims=['time'])
>>> data['time'] = np.arange(1000) / 100.0  # 100 Hz
>>> filtered = data.sps_filters.bandpass(10, 40)  # 10-40 Hz
exception mth5.timeseries.scipy_filters.DecimationWarning[source]

Bases: Warning

class mth5.timeseries.scipy_filters.FilterAccessor(darray: DataArray | Dataset)[source]

Bases: object

Accessor exposing common frequency and filtering methods.

Registered as xarray accessor under .sps_filters for both DataArray and Dataset objects.

darray[source]

The wrapped xarray object.

Type:

xarray.DataArray | xarray.Dataset

Examples

Apply filters via accessor:

>>> data.sps_filters.low(10)  # lowpass at 10 Hz
>>> data.sps_filters.bandpass(5, 15)  # bandpass 5-15 Hz
bandpass(f_low: float, f_high: float, *args: Any, **kwargs: Any) DataArray | Dataset[source]

Apply bandpass filter.

Parameters:
  • f_low (float) – Lower cutoff frequency in Hz.

  • f_high (float) – Upper cutoff frequency in Hz.

  • *args – Passed to bandpass: (order, irtype, filtfilt, apply_kwargs, in_nyq, dim).

  • **kwargs – Passed to filter design.

Returns:

Bandpass-filtered data.

Return type:

xarray.DataArray | xarray.Dataset

bandstop(f_low: float, f_high: float, *args: Any, **kwargs: Any) DataArray | Dataset[source]

Apply bandstop filter.

Parameters:
  • f_low (float) – Lower cutoff frequency in Hz.

  • f_high (float) – Upper cutoff frequency in Hz.

  • *args – Passed to bandstop: (order, irtype, filtfilt, apply_kwargs, in_nyq, dim).

  • **kwargs – Passed to filter design.

Returns:

Bandstop-filtered data.

Return type:

xarray.DataArray | xarray.Dataset

decimate(target_sample_rate: float, n_order: int = 8, dim: str | None = None) DataArray | Dataset[source]

Decimate signal.

Parameters:
  • target_sample_rate (float) – Target sample rate.

  • n_order (int, default 8) – Chebyshev filter order.

  • dim (str | None, default None) – Dimension to decimate along.

Returns:

Decimated data.

Return type:

xarray.DataArray | xarray.Dataset

detrend(trend_type: str = 'linear', dim: str | None = None) DataArray | Dataset[source]

Remove trend from data.

Parameters:
  • trend_type ({'linear', 'constant'}, default 'linear') – Type of trend to remove.

  • dim (str | None, default None) – Dimension to detrend along.

Returns:

Detrended data.

Return type:

xarray.DataArray | xarray.Dataset

property dt: float[source]

Sampling step of last axis.

property dx: ndarray[source]

Sampling steps for all axes as array.

freq(f_crit: float | list[float] | tuple[float, float], order: int | None = None, irtype: str = 'iir', filtfilt: bool = True, apply_kwargs: dict[str, Any] | None = None, in_nyq: bool = False, dim: str | None = None, **kwargs: Any) DataArray | Dataset[source]

Apply general frequency filter.

Parameters:
  • f_crit (float | list[float] | tuple[float, float]) – Critical frequency or frequencies in Hz.

  • order (int | None, default None) – Filter order.

  • irtype ({'iir', 'fir'}, default 'iir') – Impulse response type.

  • filtfilt (bool, default True) – Apply forward-backward filtering.

  • apply_kwargs (dict | None, default None) – Additional filter function kwargs.

  • in_nyq (bool, default False) – If True, f_crit is normalized by Nyquist frequency.

  • dim (str | None, default None) – Dimension along which to filter.

  • **kwargs – Passed to filter design.

Returns:

Filtered data.

Return type:

xarray.DataArray | xarray.Dataset

property fs: float[source]

Sampling frequency in inverse units of self.dt.

high(f_cutoff: float, *args: Any, **kwargs: Any) DataArray | Dataset[source]

Apply highpass filter.

Parameters:
  • f_cutoff (float) – Cutoff frequency in Hz.

  • *args – Passed to highpass: (order, irtype, filtfilt, apply_kwargs, in_nyq, dim).

  • **kwargs – Passed to filter design.

Returns:

Highpass-filtered data.

Return type:

xarray.DataArray | xarray.Dataset

low(f_cutoff: float, *args: Any, **kwargs: Any) DataArray | Dataset[source]

Apply lowpass filter.

Parameters:
  • f_cutoff (float) – Cutoff frequency in Hz.

  • *args – Passed to lowpass: (order, irtype, filtfilt, apply_kwargs, in_nyq, dim).

  • **kwargs – Passed to filter design.

Returns:

Lowpass-filtered data.

Return type:

xarray.DataArray | xarray.Dataset

resample_poly(target_sample_rate: float, pad_type: str = 'mean', dim: str | None = None) DataArray | Dataset[source]

Resample using polyphase filtering.

Parameters:
  • target_sample_rate (float) – Target sample rate.

  • pad_type (str, default 'mean') – Padding type for resampling.

  • dim (str | None, default None) – Dimension to resample along.

Returns:

Resampled data.

Return type:

xarray.DataArray | xarray.Dataset

savgol(window_length: int, polyorder: int, deriv: int = 0, delta: float | None = None, dim: str | None = None, mode: str = 'interp', cval: float = 0.0) DataArray | Dataset[source]

Apply Savitzky-Golay filter.

Parameters:
  • window_length (int) – Filter window length (positive odd integer).

  • polyorder (int) – Polynomial order (< window_length).

  • deriv (int, default 0) – Derivative order.

  • delta (float | None, default None) – Sample spacing.

  • dim (str | None, default None) – Dimension to filter along.

  • mode (str, default 'interp') – Extension mode.

  • cval (float, default 0.0) – Constant fill value.

Returns:

Filtered data.

Return type:

xarray.DataArray | xarray.Dataset

exception mth5.timeseries.scipy_filters.FilteringNaNWarning[source]

Bases: Warning

exception mth5.timeseries.scipy_filters.UnevenSamplingWarning[source]

Bases: Warning

mth5.timeseries.scipy_filters.bandpass(darray: DataArray | Dataset, f_low: float, f_high: float, *args: Any, **kwargs: Any) DataArray | Dataset[source]

Apply a bandpass filter.

Parameters:
  • darray (xarray.DataArray | xarray.Dataset) – Data to filter.

  • f_low (float) – Lower cutoff frequency in Hz.

  • f_high (float) – Upper cutoff frequency in Hz.

  • *args – Passed to frequency_filter: (order, irtype, filtfilt, apply_kwargs, in_nyq, dim).

  • **kwargs – Passed to filter design (see frequency_filter).

Returns:

Bandpass-filtered data.

Return type:

xarray.DataArray | xarray.Dataset

Examples

Keep components between 10 and 50 Hz:

>>> filtered = bandpass(data, 10, 50)
mth5.timeseries.scipy_filters.bandstop(darray: DataArray | Dataset, f_low: float, f_high: float, *args: Any, **kwargs: Any) DataArray | Dataset[source]

Apply a bandstop (notch) filter.

Parameters:
  • darray (xarray.DataArray | xarray.Dataset) – Data to filter.

  • f_low (float) – Lower cutoff frequency in Hz.

  • f_high (float) – Upper cutoff frequency in Hz.

  • *args – Passed to frequency_filter: (order, irtype, filtfilt, apply_kwargs, in_nyq, dim).

  • **kwargs – Passed to filter design (see frequency_filter).

Returns:

Bandstop-filtered data (removes frequencies between f_low and f_high).

Return type:

xarray.DataArray | xarray.Dataset

Examples

Remove 50-60 Hz powerline noise:

>>> filtered = bandstop(data, 50, 60)
mth5.timeseries.scipy_filters.decimate(darray: Dataset | DataArray, target_sample_rate: float, n_order: int = 8, dim: str | None = None) DataArray | Dataset[source]

Decimate data using Chebyshev filter and downsampling.

Applies an 8th-order Chebyshev type I filter with zero-phase filtering (sosfiltfilt) before downsampling.

Parameters:
  • darray (xarray.DataArray | xarray.Dataset) – Data to decimate.

  • target_sample_rate (float) – Target sample rate in samples per second (or per space unit).

  • n_order (int, default 8) – Order of the Chebyshev type I filter.

  • dim (str | None, default None) – Dimension to decimate along. Auto-detected for 1-D arrays.

Returns:

Decimated data with adjusted coordinates.

Return type:

xarray.DataArray | xarray.Dataset

Raises:

ValueError – If dim is None and array is not 1-D.

Warning

UserWarning

If decimation factor > 13, suggest calling decimate multiple times.

Notes

If sample_rate / target_sample_rate > 13, call decimate multiple times to avoid aliasing artifacts.

Examples

Decimate from 100 Hz to 10 Hz:

>>> decimated = decimate(data, target_sample_rate=10.0)
mth5.timeseries.scipy_filters.detrend(darray: DataArray | Dataset, dim: str | None = None, trend_type: str = 'linear') DataArray | Dataset[source]

Remove linear or constant trend from data.

Parameters:
  • darray (xarray.DataArray | xarray.Dataset) – Data to detrend.

  • dim (str | None, default None) – Dimension along which to detrend. Auto-detected for 1-D arrays.

  • trend_type ({'linear', 'constant'}, default 'linear') – Type of detrending: ‘linear’ removes linear trend via least-squares, ‘constant’ removes mean.

Returns:

Detrended data.

Return type:

xarray.DataArray | xarray.Dataset

Raises:

ValueError – If dim is None and array is not 1-D.

Examples

Remove linear trend:

>>> detrended = detrend(data, trend_type='linear')

Remove mean (DC component):

>>> demeaned = detrend(data, trend_type='constant')
mth5.timeseries.scipy_filters.frequency_filter(darray: DataArray | Dataset, f_crit: float | list[float] | tuple[float, float], order: int | None = None, irtype: str = 'iir', filtfilt: bool = True, apply_kwargs: dict[str, Any] | None = None, in_nyq: bool = False, dim: str | None = None, **kwargs: Any) DataArray | Dataset[source]

Apply a frequency filter to an xarray.

Supports IIR (infinite impulse response) and FIR (finite impulse response) filters with optional forward-backward filtering for zero-phase response.

Parameters:
  • darray (xarray.DataArray | xarray.Dataset) – Data to be filtered.

  • f_crit (float | list[float] | tuple[float, float]) – Critical frequency or frequencies (in Hz). Scalar for lowpass/highpass, 2-element for bandpass/bandstop.

  • order (int | None, default None) – Filter order. If None, uses default (8 for IIR, 29 for FIR).

  • irtype ({'iir', 'fir'}, default 'iir') – Impulse response type: ‘iir’ or ‘fir’.

  • filtfilt (bool, default True) – Apply filter forwards and backwards for zero-phase response.

  • apply_kwargs (dict | None, default None) – Additional kwargs passed to the filter function.

  • in_nyq (bool, default False) – If True, f_crit values are already normalized by Nyquist frequency.

  • dim (str | None, default None) – Dimension along which to filter. Auto-detected for 1-D arrays.

  • **kwargs – Passed to filter design function (scipy.signal.iirfilter or firwin).

Returns:

Filtered data with same structure as input.

Return type:

xarray.DataArray | xarray.Dataset

Raises:

ValueError – If irtype is not ‘iir’ or ‘fir’, or if dim is ambiguous.

Warning

FilteringNaNWarning

If input contains NaN values.

Examples

Apply a 4th-order IIR lowpass at 10 Hz:

>>> filtered = frequency_filter(data, 10, order=4, btype='low')

FIR bandpass filter from 5 to 15 Hz:

>>> filtered = frequency_filter(data, [5, 15], irtype='fir', btype='band')
mth5.timeseries.scipy_filters.get_maybe_last_dim_axis(darray: DataArray | Dataset, dim: str | None = None) tuple[str, int][source]

Get dimension name and axis index.

Parameters:
  • darray (xarray.DataArray | xarray.Dataset) – Input array.

  • dim (str | None, default None) – Dimension name. If None, uses the last dimension.

Returns:

Dimension name and corresponding axis index.

Return type:

tuple[str, int]

mth5.timeseries.scipy_filters.get_maybe_only_dim(darray: DataArray | Dataset, dim: str | None) str[source]

Determine the dimension along which to operate.

If dim is None and the array is 1-D, returns the single dimension. Otherwise, returns the provided dim.

Parameters:
  • darray (xarray.DataArray | xarray.Dataset) – An xarray DataArray or Dataset.

  • dim (str | None) – Dimension name, or None to auto-detect for 1-D arrays.

Returns:

The dimension name.

Return type:

str

Raises:

ValueError – If dim is None and the array is not 1-D.

mth5.timeseries.scipy_filters.get_sampling_step(darray: DataArray | Dataset, dim: str | None = None, rtol: float = 0.001) float[source]

Compute the sampling step along a dimension.

Automatically detects time unit (ns, us, ms, s) and scales accordingly. Issues a warning if sampling is not uniform (average vs first step mismatch).

Parameters:
  • darray (xarray.DataArray | xarray.Dataset) – Input array with coordinates.

  • dim (str | None, default None) – Dimension name. Auto-detected for 1-D arrays.

  • rtol (float, default 1e-3) – Relative tolerance for detecting uneven sampling.

Returns:

Sampling step in appropriate units (seconds for time coordinates).

Return type:

float

Raises:

ValueError – If the coordinate has fewer than 2 samples.

Warning

UnevenSamplingWarning

If average sampling step differs from first step by more than rtol.

Examples

Get sampling step from a time series:

>>> dt = get_sampling_step(data, dim='time')
>>> sample_rate = 1.0 / dt
mth5.timeseries.scipy_filters.highpass(darray: DataArray | Dataset, f_cutoff: float, *args: Any, **kwargs: Any) DataArray | Dataset[source]

Apply a highpass filter.

Parameters:
  • darray (xarray.DataArray | xarray.Dataset) – Data to filter.

  • f_cutoff (float) – Cutoff frequency in Hz.

  • *args – Passed to frequency_filter: (order, irtype, filtfilt, apply_kwargs, in_nyq, dim).

  • **kwargs – Passed to filter design (see frequency_filter).

Returns:

Highpass-filtered data.

Return type:

xarray.DataArray | xarray.Dataset

Examples

Remove components below 1 Hz:

>>> filtered = highpass(data, 1.0)
mth5.timeseries.scipy_filters.lowpass(darray: DataArray | Dataset, f_cutoff: float, *args: Any, **kwargs: Any) DataArray | Dataset[source]

Apply a lowpass filter.

Parameters:
  • darray (xarray.DataArray | xarray.Dataset) – Data to filter.

  • f_cutoff (float) – Cutoff frequency in Hz.

  • *args – Passed to frequency_filter: (order, irtype, filtfilt, apply_kwargs, in_nyq, dim).

  • **kwargs – Passed to filter design (see frequency_filter).

Returns:

Lowpass-filtered data.

Return type:

xarray.DataArray | xarray.Dataset

Examples

Remove components above 50 Hz:

>>> filtered = lowpass(data, 50)
mth5.timeseries.scipy_filters.resample_poly(darray: DataArray | Dataset, new_sample_rate: float, dim: str | None = None, pad_type: str = 'mean') DataArray | Dataset[source]

Resample using polyphase filtering.

Computes rational resampling ratio (up/down) and applies scipy.signal.resample_poly. Automatically handles coordinate updates.

Parameters:
  • darray (xarray.DataArray | xarray.Dataset) – Data to resample.

  • new_sample_rate (float) – Target sample rate.

  • dim (str | None, default None) – Dimension to resample along. Auto-detected for 1-D arrays.

  • pad_type (str, default 'mean') – Padding type passed to scipy.signal.resample_poly. Options: ‘constant’, ‘line’, ‘mean’, ‘median’, etc.

Returns:

Resampled data with updated coordinates.

Return type:

xarray.DataArray | xarray.Dataset

Raises:

ValueError – If dim is None and array is not 1-D.

Warning

UserWarning

If new sample rate is not an integer multiple of original rate.

Notes

In newer scipy versions, data is cast to float and returns float dtype.

Examples

Resample to 50 Hz:

>>> resampled = resample_poly(data, new_sample_rate=50.0)
mth5.timeseries.scipy_filters.savgol_filter(darray: DataArray | Dataset, window_length: int, polyorder: int, deriv: int = 0, delta: float | None = None, dim: str | None = None, mode: str = 'interp', cval: float = 0.0) DataArray | Dataset[source]

Apply a Savitzky-Golay filter.

Smooths data using least-squares polynomial fit over a sliding window. Can also compute derivatives.

Parameters:
  • darray (xarray.DataArray | xarray.Dataset) – Data to filter. Converted to float64 if not already float.

  • window_length (int) – Length of the filter window (number of coefficients). Must be positive odd integer.

  • polyorder (int) – Order of polynomial for fitting. Must be < window_length.

  • deriv (int, default 0) – Order of derivative to compute (0 = no differentiation).

  • delta (float | None, default None) – Sample spacing for derivative computation. Only used if deriv > 0.

  • dim (str | None, default None) – Dimension along which to filter. Auto-detected for 1-D arrays.

  • mode (str, default 'interp') – Extension mode: ‘mirror’, ‘constant’, ‘nearest’, ‘wrap’, or ‘interp’.

  • cval (float, default 0.0) – Fill value when mode=’constant’.

Returns:

Filtered data.

Return type:

xarray.DataArray | xarray.Dataset

Raises:

ValueError – If window_length is not positive odd, polyorder >= window_length, or dim is None for multi-dimensional arrays.

Examples

Smooth with 11-point window, 2nd-order polynomial:

>>> smoothed = savgol_filter(data, window_length=11, polyorder=2)

Compute first derivative:

>>> deriv1 = savgol_filter(data, 11, 2, deriv=1, delta=0.01)

mth5.timeseries.ts_filters module

time series filters

class mth5.timeseries.ts_filters.RemoveInstrumentResponse(ts, time_array, sample_interval, channel_response, **kwargs)[source]

Bases: object

Remove instrument response from the given channel response filter

The order of operations is important (if applied):

  1. detrend

  2. zero mean

  3. zero pad

  4. time window

  5. frequency window

  6. remove response

  7. undo time window

  8. bandpass

Parameters:
  • ts (np.ndarray((N,) , dtype=float)) – time series data to remove response from

  • time_array (np.ndarray((N,) , dtype=np.datetime[ns])) – time index that corresponds to the time series

  • sample_interval (float) – seconds per sample (time interval between samples)

  • channel_response – Channel response filter with all filters

included to convert from counts to physical units :type channel_response: class:mt_metadata.timeseries.filters.ChannelResponse`

kwargs

Parameters:
  • plot (boolean, default True) – to plot the calibration process [ False | True ]

  • detrend (boolean, default True) – Remove linar trend of the time series

  • zero_mean (boolean, default True) – Remove the mean of the time series

  • zero_pad (boolean, default True) – pad the time series to the next power of 2 for efficiency

  • t_window (string, default None) – Time domain window name see scipy.signal.windows for options

  • t_window_params – Time domain window parameters, parameters can be

found in scipy.signal.windows :type t_window_params: dictionary :param f_window: Frequency domain window name see scipy.signal.windows for options :type f_window: string, default None :param f_window_params: Frequency window parameters, parameters can be found in scipy.signal.windows :type f_window_params: dictionary :param bandpass: bandpass frequency and order {“low”:, “high”:, “order”:,} :type bandpass: dictionary

apply_bandpass(ts)[source]

apply a bandpass filter to the calibrated data

Parameters:

ts (np.ndarray) – calibrated time series

Returns:

bandpassed time series

Return type:

np.ndarray

apply_detrend(ts)[source]

Detrend time series using scipy.detrend(‘linear’)

Parameters:

ts (np.ndarray) – input time series

Returns:

detrended time series

Return type:

np.ndarray

apply_f_window(data)[source]

Apply a frequency domain window. Get the available windows from scipy.signal.windows

Need to create a window twice the size of the input because we are only taking the rfft which gives just half the spectra and then take only half the window

Parameters:

data (np.ndarray) – input spectra

Returns:

windowed spectra

Return type:

np.ndarray

apply_t_window(ts)[source]

Apply a window in the time domain. Get the available windows from scipy.signal.windows

Parameters:

ts (np.ndarray) – input time series

Returns:

windowed time series

Return type:

np.ndarray

apply_zero_mean(ts)[source]

Remove the mean from the time series

Parameters:

ts (np.ndarray) – input time series

Returns:

zero mean time series

Return type:

np.ndarray

apply_zero_pad(ts)[source]

zero pad to power of 2, at the end of the time series to make the FFT more efficient

Parameters:

ts (np.ndarray) – input time series

Returns:

zero padded time series

Return type:

np.ndarray

static get_window(window, window_params, size)[source]

Get window from scipy.signal

Parameters:
  • window (string) – name of the window

  • window_params (dictionary) – dictionary of window parameters

  • size (integer) – number of points in the window

Returns:

window function

Return type:

class:scipy.signal

remove_instrument_response(operation='divide', include_decimation=None, include_delay=None, filters_to_remove=[])[source]

Remove instrument response following the recipe provided

Returns:

calibrated time series

Return type:

np.ndarray

mth5.timeseries.ts_filters.adaptive_notch_filter(bx, df=100, notches=[50, 100], notchradius=0.5, freqrad=0.9, rp=0.1, dbstop_limit=5.0)[source]
Parameters:
  • bx (np.ndarray) – time series to filter

  • df (float, optional) – sample rate in samples per second, defaults to 100

  • notches (list, optional) – list of frequencies to locate notches at in Hz, defaults to [50, 100]

  • notchradius (float, optional) – notch radius, defaults to 0.5

  • freqrad (float, optional) – radius to search for a peak at the notch frequency, defaults to 0.9

  • rp (float, optional) – ripple of Chebyshev type 1 filter, lower numbers means less ripples, defaults to 0.1

  • dbstop_limit (float, optional) – limits the difference between the peak at the notch and surrounding spectra. Any difference above dbstop_limit will be filtered, anything less will not, defaults to 5.0

Returns:

notch filtered data

Return type:

np.ndarray

Returns:

list of notch frequencies

Return type:

list

Example

>>> import RemovePeriodicNoise_Kate as rmp
>>> # make a variable for the file to load in
>>> fn = r"/home/MT/mt01_20130101_000000.BX"
>>> # load in file, if the time series is not an ascii file
>>> # might need to add keywords to np.loadtxt or use another
>>> # method to read in the file
>>> bx = np.loadtxt(fn)
>>> # create a list of frequencies to filter out
>>> freq_notches = [50, 150, 200]
>>> # filter data
>>> bx_filt, filt_lst = rmp.adaptiveNotchFilter(bx, df=100.
>>> ...                                         notches=freq_notches)
>>> #save the filtered data into a file
>>> np.savetxt(r"/home/MT/Filtered/mt01_20130101_000000.BX", bx_filt)

Notes

Most of the time the default parameters work well, the only thing you need to change is the notches and perhaps the radius. I would test it out with a few time series to find the optimum parameters. Then make a loop over all you time series data. Something like

>>> import os
>>> dirpath = r"/home/MT"
>>> #make a director to save filtered time series
>>> save_path = r"/home/MT/Filtered"
>>> if not os.path.exists(save_path):
>>>     os.mkdir(save_path)
>>> for fn in os.listdir(dirpath):
>>>     bx = np.loadtxt(os.path.join(dirpath, fn)
>>>     bx_filt, filt_lst = rmp.adaptiveNotchFilter(bx, df=100.
>>>     ...                                         notches=freq_notches)
>>>     np.savetxt(os.path.join(save_path, fn), bx_filt)
mth5.timeseries.ts_filters.butter_bandpass(lowcut, highcut, sample_rate, order=5)[source]

Butterworth bandpass filter using scipy.signal

Transforms band corners to angular frequencies

Parameters:
  • lowcut (float) – low cut frequency in Hz (3dB point)

  • highcut (float) – high cut frequency in Hz (3dB point)

  • sample_rate (float) – Sample rate

  • order (int, optional) – Butterworth order, defaults to 5

Returns:

SOS scipy.signal format

Return type:

scipy.signal.SOS?

mth5.timeseries.ts_filters.butter_bandpass_filter(data, lowcut, highcut, fs, order=5)[source]
Parameters:
  • data (np.ndarray) – 1D time series data

  • lowcut (float) – low cut frequency in Hz

  • highcut (float) – high cut frequency in Hz

  • fs (float) – Sample rate

  • order (int, optional) – Butterworth order, defaults to 5

Returns:

filtered data

Return type:

np.ndarray

mth5.timeseries.ts_filters.low_pass(data, low_pass_freq, cutoff_freq, sample_rate)[source]
Parameters:
  • data (np.ndarray) – 1D time series data

  • low_pass_freq (float) – low pass frequency in Hz

  • cutoff_freq (float) – cut off frequency in Hz

  • sampling_rate (float) – Sample rate in samples per second

Returns:

lowpass filtered data

Return type:

np.ndarray

mth5.timeseries.ts_filters.zero_pad(input_array, power=2, pad_fill=0)[source]
Parameters:
  • input_array (np.ndarray) – 1D array

  • power – base power to used to pad to, defaults to 2 which is optimal

for the FFT :type power: int, optional :param pad_fill: fill value for padded values, defaults to 0 :type pad_fill: float, optional :return: zero padded array :rtype: np.ndarray

mth5.timeseries.ts_helpers module

Time series helper functions.

Provides utilities for time series processing including decimation planning, significant figure counting, and datetime coordinate generation.

Notes

Created on Wed Mar 29 14:27:22 2023

Author: jpeacock

mth5.timeseries.ts_helpers.get_decimation_sample_rates(old_sample_rate: float, new_sample_rate: float, max_decimation: int) list[float][source]

Compute decimation steps from old to new sample rate.

Generates a list of intermediate sample rates to decimate from old_sample_rate to new_sample_rate without exceeding max_decimation at any single step. Uses geometric series to distribute decimation evenly.

Parameters:
  • old_sample_rate (float) – Original sample rate in samples per second.

  • new_sample_rate (float) – Target sample rate in samples per second.

  • max_decimation (int) – Maximum allowed decimation factor per step.

Returns:

List of intermediate sample rates for multi-stage decimation. If total decimation <= max_decimation, returns [new_sample_rate].

Return type:

list[float]

Examples

Single-step decimation (factor <= max):

>>> get_decimation_sample_rates(100, 10, 13)
[10]

Multi-step decimation (factor > max):

>>> get_decimation_sample_rates(1000, 10, 13)
[77, 10]
>>> # Decimates 1000 -> 77 -> 10

Three-step decimation:

>>> get_decimation_sample_rates(10000, 10, 13)
[770, 60, 10]
mth5.timeseries.ts_helpers.make_dt_coordinates(start_time: str | MTime | None, sample_rate: float | None, n_samples: int, end_time: str | MTime | None = None) DatetimeIndex[source]

Generate datetime coordinates for time series data.

Creates a pandas DatetimeIndex with proper rounding based on sample rate precision. Handles edge cases like invalid sample rates and missing start times with sensible defaults.

Parameters:
  • start_time (str | MTime | None) – Start time in ISO format or MTime object. If None, defaults to ‘1980-01-01T00:00:00’.

  • sample_rate (float | None) – Sample rate in samples per second. If None or 0, defaults to 1 and issues a warning.

  • n_samples (int) – Number of samples in the time series. Must be >= 1.

  • end_time (str | MTime | None, default None) – End time in ISO format or MTime object. If None, computed from start_time, sample_rate, and n_samples.

Returns:

DatetimeIndex with n_samples timestamps, rounded to appropriate precision (ms, us, or ns) based on significant figures.

Return type:

pandas.DatetimeIndex

Raises:

ValueError – If n_samples < 1.

Warning

UserWarning

If sample_rate is invalid (None or 0) or start_time is None.

Notes

  • End time is calculated as: start_time + (n_samples - 1) / sample_rate

  • Rounding precision determined by significant figures in start time and sample period (1 / sample_rate): * < 3 sig figs: millisecond rounding * 3-5 sig figs: microsecond rounding * 6-8 sig figs: nanosecond rounding * >= 9 sig figs: no rounding

Examples

Create 100 Hz time index for 1000 samples:

>>> dt = make_dt_coordinates('2023-01-01T00:00:00', 100.0, 1000)
>>> len(dt)
1000
>>> dt[0]
Timestamp('2023-01-01 00:00:00')
>>> dt[-1]
Timestamp('2023-01-01 00:00:09.990000')

With explicit end time:

>>> dt = make_dt_coordinates(
...     '2023-01-01T00:00:00',
...     10.0,
...     100,
...     end_time='2023-01-01T00:00:09.9'
... )

Handle invalid sample rate (uses default of 1 Hz):

>>> dt = make_dt_coordinates('2023-01-01T00:00:00', None, 10)
# Warning: Need to input a valid sample rate...

mth5.timeseries.xarray_helpers module

Module containing helper functions for working with xarray objects.

mth5.timeseries.xarray_helpers.covariance_xr(X: DataArray, aweights: ndarray | None = None, bias: bool | None = False, rowvar: bool | None = False) DataArray[source]

Compute the covariance matrix with numpy.cov.

Parameters:
  • X (xarray.core.dataarray.DataArray) – Multivariate time series as an xarray

  • aweights (array_like, optional) – Passthrough param for np.cov. 1-D array of observation vector weights. These relative weights are typically large for observations considered “important” and smaller for observations considered less “important”. If ddof=0 the array of weights can be used to assign probabilities to observation vectors.

  • bias (bool) – Passthrough param for np.cov. Default normalization (False) is by (N - 1), where N is the number of observations given (unbiased estimate). If bias is True, then normalization is by N. These values can be overridden by using the keyword ddof in numpy versions >= 1.5.

  • rowvar (bool) – Passthrough param for np.cov. If rowvar is True (default), then each row represents a variable, with observations in the columns. Otherwise, the relationship is transposed: each column represents a variable, while the rows contain observations.

Returns:

  • S (xarray.DataArray) – The covariance matrix of the data in xarray form.

  • Development Notes – In case of ValueError: conflicting sizes for dimension ‘channel_1’, this likely means the bool for rowvar should be flipped.

mth5.timeseries.xarray_helpers.initialize_xrda_1d(channels: list, dtype: type | None = None, value: complex | float | bool | None = 0.0) DataArray[source]

Returns a 1D xr.DataArray with variable “channel”, having values channels named by the input list.

Parameters:
  • channels (list) – The channels in the multivariate array

  • dtype (type, optional) – The datatype to initialize the array. Common cases are complex, float, and bool

  • value (Union[complex, float, bool], optional) – The default value to assign the array

Returns:

xrda – An xarray container for the channels, initialized to zeros.

Return type:

xarray.core.dataarray.DataArray

mth5.timeseries.xarray_helpers.initialize_xrda_2d(variables: list, coords: dict, dtype: type = <class 'complex'>, value: int | float = 0.0) Dataset[source]

Initialize a 3D xarray DataArray with dimensions from coords plus ‘variable’.

Parameters:
  • variables (list) – List of variable names for the additional dimension.

  • coords (dict) – Dictionary of coordinates for the dataset dimensions.

  • dtype (type, optional) – Data type for the array, by default complex.

  • value (int or float, optional) – Value to initialize the array with, by default 0.

Returns:

A 3D DataArray with dimensions from coords plus ‘variable’.

Return type:

xr.DataArray

mth5.timeseries.xarray_helpers.initialize_xrda_2d_cov(channels, dtype=<class 'complex'>, value: complex | float | bool | None = 0)[source]

TODO: consider changing nomenclature from dims=[“channel_1”, “channel_2”], to dims=[“variable_1”, “variable_2”], to be consistent with initialize_xrda_1d

Parameters
channels: list

The channels in the multivariate array. The covariance matrix will be square with dimensions len(channels) x len(channels).

dtype: type

The datatype to initialize the array. Common cases are complex, float, and bool

value: Union[complex, float, bool]

The default value to assign the array

Returns
xrda: xarray.core.dataarray.DataArray

An xarray container for the channel covariances, initialized to zeros. The matrix is square with dimensions len(channels) x len(channels).

mth5.timeseries.xarray_helpers.initialize_xrds_2d(variables: list, coords: dict, dtype: type | None = <class 'complex'>, value: complex | float | bool | None = 0) Dataset[source]

Returns a 2D xr.Dataset with the given variables and coordinates.

Parameters:
  • variables (list) – List of variable names to create in the dataset

  • coords (dict) – Dictionary of coordinates for the dataset dimensions

  • dtype (type, optional) – The datatype to initialize the arrays. Common cases are complex, float, and bool

  • value (Union[complex, float, bool], optional) – The default value to assign the arrays

Returns:

xrds – A 2D xarray Dataset with dimensions from coords

Return type:

xr.Dataset

Module contents

class mth5.timeseries.ChannelTS(channel_type: str = 'auxiliary', data: ndarray | DataFrame | Series | DataArray | list | tuple | None = None, channel_metadata: Electric | Magnetic | Auxiliary | dict | None = None, station_metadata: Station | dict | None = None, run_metadata: Run | dict | None = None, survey_metadata: Survey | dict | None = None, **kwargs: Any)[source]

Bases: object

Time series container for a single MT channel with full metadata.

Stores equally-spaced time series data in an xarray.DataArray with a time coordinate index. Integrates comprehensive metadata from Survey/Station/Run/Channel hierarchy and supports calibration, resampling, merging, and format conversions.

Parameters:
  • channel_type ({'electric', 'magnetic', 'auxiliary'}, default 'auxiliary') – Type of the channel.

  • data (array-like, optional) – Time series data (numpy array, pandas DataFrame/Series, xarray.DataArray).

  • channel_metadata (mt_metadata.timeseries.Electric | Magnetic | Auxiliary | dict, optional) – Channel-specific metadata.

  • station_metadata (mt_metadata.timeseries.Station | dict, optional) – Station metadata.

  • run_metadata (mt_metadata.timeseries.Run | dict, optional) – Run metadata.

  • survey_metadata (mt_metadata.timeseries.Survey | dict, optional) – Survey metadata.

  • **kwargs – Additional attributes to set on the object.

ts

The time series data array.

Type:

numpy.ndarray

sample_rate

Sample rate in samples per second.

Type:

float

start

Start time (UTC).

Type:

MTime

end

End time (UTC), derived from start + duration.

Type:

MTime

n_samples

Number of samples.

Type:

int

component

Component name (e.g., ‘ex’, ‘hy’, ‘temperature’).

Type:

str

channel_response

Full instrument response filter chain.

Type:

ChannelResponse

Notes

  • End time is a derived property and cannot be set directly.

  • Leverages xarray for efficient interpolation, resampling, and groupby operations.

  • Metadata follows mt_metadata standards with automatic time period updates.

Examples

Create an auxiliary channel with synthetic data:

>>> from mth5.timeseries import ChannelTS
>>> import numpy as np
>>> ts_obj = ChannelTS('auxiliary')
>>> ts_obj.sample_rate = 8
>>> ts_obj.start = '2020-01-01T12:00:00+00:00'
>>> ts_obj.ts = np.random.randn(4096)
>>> ts_obj.station_metadata.id = 'MT001'
>>> ts_obj.run_metadata.id = 'MT001a'
>>> ts_obj.component = 'temperature'
>>> print(ts_obj)

Calibrate and remove instrument response:

>>> calibrated = ts_obj.remove_instrument_response()
>>> calibrated.channel_metadata.units
property channel_metadata: Electric | Magnetic | Auxiliary

Channel metadata.

Returns:

Channel metadata from the first channel in the run.

Return type:

mt_metadata.timeseries.Electric | Magnetic | Auxiliary

property channel_response

Full channel response filter

Returns:

full channel response filter

Return type:

mt_metadata.timeseries.filters.ChannelResponse

property channel_type: str

Channel type.

Returns:

Channel type: ‘Electric’, ‘Magnetic’, or ‘Auxiliary’.

Return type:

str

property component
compute_sample_rate()[source]

Two cases, high_frequency (HF) data and not HF data.

# Original comment about the HF case: Taking the median(diff(timestamps)) is more accurate for high sample rates, the way pandas.date_range rounds nanoseconds is not consistent between samples, therefore taking the median provides better results if the time series is long this can be inefficient so test first

copy(data: bool = True) ChannelTS[source]

Create a copy of the ChannelTS object.

Parameters:

data (bool, default True) – Include data in the copy (True) or only metadata (False).

Returns:

Copy of the channel.

Return type:

ChannelTS

Examples

Copy metadata structure without data:

>>> ch_copy = ts_obj.copy(data=False)
decimate(new_sample_rate, inplace=False, max_decimation=8)[source]

decimate the data by using scipy.signal.decimate

Parameters:

dec_factor (int) – decimation factor

  • refills ts.data with decimated data and replaces sample_rate

property end

MTime object

from_obspy_trace(obspy_trace)[source]

Fill data from an obspy.core.Trace

Parameters:

obspy_trace (obspy.core.trace) – Obspy trace object

get_calibrated_units()[source]

Follows the FDSN standard which has the filter stages starting with physical units to digital counts.

The channel_response is expected to have a list of filter “stages” of which the first stage has input units corresponding to the the physical quantity that the instrument measures, and the last is normally counts.

channel_response can be viewed as the chaining together of all of these filters.

Thus it is normal for channel_response.units_out will be in the same units as the archived raw time series, and for the units after the response is corrected for will be the units_in of

The units of the channel metadata are compared to the input and output units of the channel_response.

Returns:

tuple, calibration_operation, either “mulitply” or divide”, and a string for calibrated units

Return type:

tuple (of two strings_

get_calibration_operation()[source]
get_sample_rate_supplied_at_init(channel_metadata: Electric | Magnetic | Auxiliary | dict | None) float | None[source]

Extract sample_rate from channel_metadata if available.

Parameters:

channel_metadata (mt_metadata.timeseries.Electric | Magnetic | Auxiliary | dict | None) – Metadata that may contain a sample_rate field.

Returns:

Sample rate if found, otherwise None.

Return type:

float | None

Notes

Supports nested dict structures like {"electric": {"sample_rate": 8.0}}.

get_slice(start, end=None, n_samples=None)[source]

Get a slice from the time series given a start and end time.

Looks for >= start & <= end

Uses loc to be exact with milliseconds

Parameters:
  • start (string, MTime) – start time of the slice

  • end (string, MTime) – end time of the slice

  • n_samples (integer) – number of sample to get after start time

Returns:

slice of the channel requested

Return type:

ChannelTS

has_data()[source]

check to see if there is an index in the time series

is_high_frequency(threshold_dt=0.0001)[source]

Quasi hard-coded condition to check if data are logged at more than 10kHz can be parameterized in future

merge(other, gap_method='slinear', new_sample_rate=None, resample_method='poly')[source]

merg two channels or list of channels together in the following steps

  1. xr.combine_by_coords([original, other])

  2. compute monotonic time index

  3. reindex(new_time_index, method=gap_method)

If you want a different method or more control use merge

Parameters:

other (mth5.timeseries.ChannelTS) – Another channel

Raises:
  • TypeError – If input is not a ChannelTS

  • ValueError – if the components are different

Returns:

Combined channel with monotonic time index and same metadata

Return type:

mth5.timeseries.ChannelTS

property n_samples

number of samples

plot()[source]

Simple plot of the data

Returns:

figure object

Return type:

matplotlib.figure

plot_spectra(spectra_type='welch', window_length=4096, **kwargs)[source]
Parameters:
  • spectra_type (string, optional) – spectra type, defaults to “welch”

  • window_length (int, optional) – window length of the welch method should be a power of 2, defaults to 2 ** 12

  • **kwargs

    DESCRIPTION

remove_instrument_response(include_decimation=False, include_delay=False, **kwargs)[source]

Remove instrument response from the given channel response filter

The order of operations is important (if applied):

  1. detrend

  2. zero mean

  3. zero pad

  4. time window

  5. frequency window

  6. remove response

  7. undo time window

  8. bandpass

Parameters:
  • include_decimation (bool, optional) – Include decimation in response, defaults to True

  • include_delay (bool, optional) – include delay in complex response, defaults to False

kwargs

Parameters:
  • plot (boolean, default True) – to plot the calibration process [ False | True ]

  • detrend (boolean, default True) – Remove linar trend of the time series

  • zero_mean (boolean, default True) – Remove the mean of the time series

  • zero_pad (boolean, default True) – pad the time series to the next power of 2 for efficiency

  • t_window (string, default None) – Time domain windown name see scipy.signal.windows for options

  • t_window_params – Time domain window parameters, parameters can be

found in scipy.signal.windows :type t_window_params: dictionary :param f_window: Frequency domain windown name see scipy.signal.windows for options :type f_window: string, defualt None :param f_window_params: Frequency window parameters, parameters can be found in scipy.signal.windows :type f_window_params: dictionary :param bandpass: bandpass freequency and order {“low”:, “high”:, “order”:,} :type bandpass: dictionary

resample_poly(new_sample_rate, pad_type='mean', inplace=False)[source]

Use scipy.signal.resample_poly to resample data while using an FIR filter to remove aliasing.

Parameters:
  • new_sample_rate (TYPE) – DESCRIPTION

  • pad_type (TYPE, optional) – DESCRIPTION, defaults to “mean”

Returns:

DESCRIPTION

Return type:

TYPE

property run_metadata: Run

Run metadata.

Returns:

Run metadata from the first run in the station.

Return type:

mt_metadata.timeseries.Run

property sample_interval

Sample interval = 1 / sample_rate

Returns:

sample interval as time distance between time samples

Return type:

float

property sample_rate

sample rate in samples/second

property start

MTime object

property station_metadata: Station

Station metadata.

Returns:

Station metadata from the first station in the survey.

Return type:

mt_metadata.timeseries.Station

property survey_metadata: Survey

Survey metadata.

Returns:

Survey metadata with updated keys.

Return type:

mt_metadata.timeseries.Survey

property time_index: ndarray

Time index as a numpy array.

Returns:

Array of datetime64[ns] timestamps.

Return type:

numpy.ndarray

to_obspy_trace(network_code=None, encoding=None)[source]

Convert the time series to an obspy.core.trace.Trace object. This will be helpful for converting between data pulled from IRIS and data going into IRIS.

Parameters:

network_code (string) – two letter code provided by FDSN DMC

Returns:

DESCRIPTION

Return type:

TYPE

to_xarray()[source]

Returns a xarray.DataArray object of the channel timeseries this way metadata from the metadata class is updated upon return.

Returns:

Returns a xarray.DataArray object of the channel timeseries

this way metadata from the metadata class is updated upon return. :rtype: xarray.DataArray

>>> import numpy as np
>>> from mth5.timeseries import ChannelTS
>>> ex = ChannelTS("electric")
>>> ex.start = "2020-01-01T12:00:00"
>>> ex.sample_rate = 16
>>> ex.ts = np.random.rand(4096)
property ts: ndarray

Time series data as a numpy array.

Returns:

The time series data.

Return type:

numpy.ndarray

welch_spectra(window_length=4096, **kwargs)[source]

get welch spectra

Parameters:
  • window_length (TYPE) – DESCRIPTION

  • **kwargs

    DESCRIPTION

Returns:

DESCRIPTION

Return type:

TYPE

class mth5.timeseries.RunTS(array_list: list[ChannelTS] | list[DataArray] | Dataset | None = None, run_metadata: Run | dict | None = None, station_metadata: Station | dict | None = None, survey_metadata: Survey | dict | None = None)[source]

Bases: object

Container for MT time series data from a single run.

Holds all run time series in one aligned xarray Dataset with channels as data variables and time as the coordinate. Manages metadata for survey, station, and run levels.

Parameters:
  • array_list (list[ChannelTS] | list[xr.DataArray] | xr.Dataset | None, optional) – List of ChannelTS objects, xarray DataArrays, or an xarray Dataset containing the time series data. All channels will be aligned to a common time index.

  • run_metadata (timeseries.Run | dict | None, optional) – Metadata for the run. Can be a Run object or dictionary.

  • station_metadata (timeseries.Station | dict | None, optional) – Metadata for the station. Can be a Station object or dictionary.

  • survey_metadata (timeseries.Survey | dict | None, optional) – Metadata for the survey. Can be a Survey object or dictionary.

dataset

xarray Dataset containing all channel data with time coordinate

Type:

xr.Dataset

survey_metadata

Survey-level metadata

Type:

timeseries.Survey

station_metadata

Station-level metadata

Type:

timeseries.Station

run_metadata

Run-level metadata

Type:

timeseries.Run

filters

Dictionary of channel response filters keyed by filter name

Type:

dict[str, Filter]

sample_rate

Sample rate in samples per second

Type:

float

channels

List of channel names in the dataset

Type:

list[str]

Examples

Create an empty RunTS:

>>> from mth5.timeseries import RunTS
>>> run = RunTS()

Create RunTS from ChannelTS objects:

>>> from mth5.timeseries import ChannelTS, RunTS
>>> ex = ChannelTS('electric', data=ex_data,
...                channel_metadata={'component': 'ex'})
>>> ey = ChannelTS('electric', data=ey_data,
...                channel_metadata={'component': 'ey'})
>>> run = RunTS(array_list=[ex, ey])
>>> print(run.channels)
['ex', 'ey']

Access individual channels:

>>> ex_channel = run.ex  # Returns ChannelTS object
>>> print(ex_channel.sample_rate)
256.0

See also

ChannelTS

Individual channel time series container

Notes

When multiple channels are provided with different start/end times, they will be automatically aligned using the earliest start and latest end times, with NaN values filling gaps.

add_channel(channel: DataArray | ChannelTS) None[source]

Add a channel to the dataset.

The channel must have the same sample rate and time coordinates that are compatible with the existing dataset. If start times don’t match, NaN values will be placed where timing doesn’t align.

Parameters:

channel (xr.DataArray | ChannelTS) – A channel as an xarray DataArray or ChannelTS object to add.

Raises:

ValueError – If the channel has a different sample rate than the run, or if the input is not a DataArray or ChannelTS.

Examples

Add a ChannelTS:

>>> hz = ChannelTS('magnetic', data=hz_data,
...                channel_metadata={'component': 'hz'})
>>> run.add_channel(hz)
>>> print(run.channels)
['ex', 'ey', 'hx', 'hy', 'hz']

Add an xarray DataArray:

>>> import xarray as xr
>>> data_array = xr.DataArray(data, coords={'time': times})
>>> run.add_channel(data_array)
calibrate(**kwargs) RunTS[source]

Remove instrument response from all channels.

Applies the channel response filters to calibrate each channel, creating a new run with calibrated data.

Parameters:

**kwargs – Additional keyword arguments passed to each channel’s remove_instrument_response method.

Returns:

New RunTS object with calibrated channels.

Return type:

RunTS

Examples

>>> calibrated_run = run.calibrate()
>>> # Calibration typically converts from counts to physical units

See also

ChannelTS.remove_instrument_response

Calibrate single channel

property channels: list[str]

List of channel names in the dataset.

Returns:

List of channel component names (e.g., [‘ex’, ‘ey’, ‘hx’]).

Return type:

list[str]

Examples

>>> print(run.channels)
['ex', 'ey', 'hx', 'hy', 'hz']
copy(data: bool = True) RunTS[source]

Create a copy of the RunTS object.

Parameters:

data (bool, optional) – If True, copy the data along with timeseries. If False, only copy the metadata (default is True).

Returns:

A copy of the RunTS object.

Return type:

RunTS

Examples

Create a copy with data:

>>> run_copy = run.copy()

Create a metadata-only copy:

>>> run_meta = run.copy(data=False)
>>> print(run_meta.has_data())
False
property dataset: Dataset

The xarray Dataset containing all channel data.

Returns:

Dataset with channels as data variables and time as coordinate.

Return type:

xr.Dataset

Examples

>>> print(run.dataset)
<xarray.Dataset>
Dimensions:  (time: 4096)
Coordinates:
  * time     (time) datetime64[ns] ...
Data variables:
    ex       (time) float64 ...
    ey       (time) float64 ...
decimate(new_sample_rate: float, inplace: bool = False, max_decimation: int = 8) RunTS | None[source]

Decimate data to a new sample rate using multi-stage decimation.

Applies FIR filtering and downsampling in multiple stages to achieve the target sample rate while preventing aliasing.

Parameters:
  • new_sample_rate (float) – Target sample rate in samples per second.

  • inplace (bool, optional) – If True, modify the current run. If False, return a new run (default is False).

  • max_decimation (int, optional) – Maximum decimation factor for each stage (default is 8).

Returns:

If inplace is False, returns new decimated RunTS. Otherwise None.

Return type:

RunTS | None

Examples

Decimate from 256 Hz to 1 Hz:

>>> decimated_run = run.decimate(1.0)
>>> print(decimated_run.sample_rate)
1.0

Decimate in place:

>>> run.decimate(16.0, inplace=True)
>>> print(run.sample_rate)
16.0

Notes

NaN values are filled with 0 before decimation to prevent NaN propagation. Multi-stage decimation is used to maintain signal quality and prevent aliasing.

See also

resample_poly

Alternative resampling using polyphase filtering

resample

Simple resampling without anti-aliasing

property end: MTime

End time of the run in UTC.

Returns:

End time from the dataset if data exists, otherwise from run_metadata.

Return type:

MTime

Examples

>>> print(run.end)
2020-01-01T01:00:00+00:00
property filters: dict[str, ChannelResponse]

Dictionary of channel response filters.

Returns:

Dictionary keyed by filter name containing ChannelResponse objects.

Return type:

dict[str, ChannelResponse]

Examples

>>> print(run.filters.keys())
dict_keys(['v_to_counts', 'dipole_100m'])
from_obspy_stream(obspy_stream: Stream, run_metadata: Run | None = None) None[source]

Get a run from an obspy.core.stream which is a list of obspy.core.Trace objects.

Parameters:

obspy_stream (obspy.core.Stream) – Obspy Stream object

Development Notes:
  • There is a baked in assumption here that the channel nomenclature in obspy is e1,e2,h1,h2,h3 and we want to convert to mth5 conventions ex,ey,hx,hy,hz. This should be made more flexible in the future.

  • A bug was found that was creating channels e1, ex, ey in the same run when reading from obspy – this is fixed here by renaming the components and a workaround to reset the station’s channels_recorded list.

get_slice(start: str | MTime, end: str | MTime | None = None, n_samples: int | None = None) RunTS[source]

Extract a time slice from the run.

Gets a chunk of data from the run, finding the closest points to the given parameters. Uses pandas slice_indexer for robust slicing.

Parameters:
  • start (str | MTime) – Start time of the slice (ISO format string or MTime object).

  • end (str | MTime | None, optional) – End time of the slice. Required if n_samples not provided.

  • n_samples (int | None, optional) – Number of samples to get. Required if end not provided.

Returns:

New RunTS object containing the requested slice with copies of metadata and filters.

Return type:

RunTS

Raises:

ValueError – If neither end nor n_samples is provided.

Examples

Get slice by start and end times:

>>> slice1 = run.get_slice('2020-01-01T00:00:00',
...                         '2020-01-01T00:01:00')
>>> print(slice1.start, slice1.end)

Get slice by start time and number of samples:

>>> slice2 = run.get_slice('2020-01-01T00:00:00', n_samples=1024)
>>> print(len(slice2.dataset.time))
1024

Notes

Uses pandas slice_indexer which handles near-matches better than xarray’s native slicing. The actual slice may be slightly adjusted to match available data points.

has_data() bool[source]

Check if the RunTS contains any data.

Returns:

True if channels with data exist, False otherwise.

Return type:

bool

Examples

>>> run = RunTS()
>>> print(run.has_data())
False
>>> run.add_channel(ex_channel)
>>> print(run.has_data())
True
merge(other: RunTS | list[RunTS], gap_method: str = 'slinear', new_sample_rate: float | None = None, resample_method: str = 'poly') RunTS[source]

Merge multiple runs into a single run.

Combines this run with one or more other runs, optionally resampling to a common sample rate and filling gaps with interpolation.

Parameters:
  • other (RunTS | list[RunTS]) – Another RunTS object or list of RunTS objects to merge.

  • gap_method (str, optional) – Interpolation method for filling gaps: ‘linear’, ‘nearest’, ‘zero’, ‘slinear’, ‘quadratic’, ‘cubic’ (default is ‘slinear’).

  • new_sample_rate (float | None, optional) – If provided, all runs will be resampled to this rate before merging. If None, uses the sample rate of the first run.

  • resample_method (str, optional) – Resampling method if new_sample_rate is provided: ‘decimate’ or ‘poly’ (default is ‘poly’).

Returns:

New merged RunTS object with monotonic time index.

Return type:

RunTS

Raises:

TypeError – If other is not a RunTS or list of RunTS objects.

Examples

Merge two runs:

>>> run1 = RunTS(array_list=[ex1, ey1])
>>> run2 = RunTS(array_list=[ex2, ey2])
>>> merged = run1.merge(run2)

Merge multiple runs with resampling:

>>> runs = [run1, run2, run3]
>>> merged = run1.merge(runs, new_sample_rate=1.0,
...                     resample_method='poly')

Notes

The merge process:

  1. Optionally resample all runs to common sample rate

  2. Combine datasets using xr.combine_by_coords

  3. Create monotonic time index spanning full range

  4. Interpolate to new index filling gaps

  5. Merge all filter dictionaries

Metadata is taken from the first run (self).

See also

__add__

Simple merging with + operator

plot(color_map: dict[str, tuple[float, float, float]] | None = None, channel_order: list[str] | None = None) Figure[source]

Plot all channels as time series.

Creates a multi-panel figure with each channel in its own subplot, sharing a common time axis.

Parameters:
  • color_map (dict[str, tuple[float, float, float]] | None, optional) –

    Dictionary mapping channel names to RGB color tuples (values 0-1). Default colors:

    • ex: (1, 0.2, 0.2) - red

    • ey: (1, 0.5, 0) - orange

    • hx: (0, 0.5, 1) - blue

    • hy: (0.5, 0.2, 1) - purple

    • hz: (0.2, 1, 1) - cyan

  • channel_order (list[str] | None, optional) – Order of channels from top to bottom. If None, uses order from self.channels.

Returns:

Figure object containing the plot.

Return type:

matplotlib.figure.Figure

Examples

Plot with default settings:

>>> fig = run.plot()
>>> fig.savefig('timeseries.png')

Plot with custom colors and order:

>>> colors = {'ex': (1, 0, 0), 'ey': (0, 1, 0)}
>>> fig = run.plot(color_map=colors, channel_order=['ey', 'ex'])

Warning

May be slow for large datasets (millions of points). Consider using get_slice() first to plot a subset.

plot_spectra(spectra_type: str = 'welch', color_map: dict[str, tuple[float, float, float]] | None = None, **kwargs) Figure[source]

Plot power spectral density for all channels.

Computes and plots the power spectrum of each channel on a single log-log plot with period on x-axis.

Parameters:
  • spectra_type (str, optional) – Type of spectral estimate to compute. Currently only ‘welch’ is supported (default is ‘welch’).

  • color_map (dict[str, tuple[float, float, float]] | None, optional) – Dictionary mapping channel names to RGB color tuples (values 0-1). Uses same defaults as plot().

  • **kwargs – Additional keyword arguments passed to the spectra computation method (e.g., nperseg, window for Welch’s method).

Returns:

Figure object containing the spectra plot.

Return type:

matplotlib.figure.Figure

Examples

Plot spectra with default settings:

>>> fig = run.plot_spectra()

Plot with custom Welch parameters:

>>> fig = run.plot_spectra(nperseg=1024, window='hann')

Notes

The plot shows:

  • Period (seconds) on bottom x-axis

  • Frequency (Hz) on top x-axis

  • Power (dB) on y-axis

See also

ChannelTS.welch_spectra

Compute Welch power spectrum

resample(new_sample_rate: float, inplace: bool = False) RunTS | None[source]

Resample data to a new sample rate using nearest-neighbor method.

Simple resampling without anti-aliasing filtering. Use decimate or resample_poly for better quality when downsampling.

Parameters:
  • new_sample_rate (float) – Target sample rate in samples per second.

  • inplace (bool, optional) – If True, modify current run. If False, return new run (default is False).

Returns:

If inplace is False, returns new resampled RunTS. Otherwise None.

Return type:

RunTS | None

Examples

>>> resampled_run = run.resample(128.0)
>>> print(resampled_run.sample_rate)
128.0

Warning

This method does not apply anti-aliasing filtering. For downsampling, consider using decimate() or resample_poly() instead.

See also

decimate

Proper downsampling with anti-aliasing

resample_poly

High-quality resampling with polyphase filtering

resample_poly(new_sample_rate: float, pad_type: str = 'mean', inplace: bool = False) RunTS | None[source]

Resample data using polyphase filtering.

Uses scipy.signal.resample_poly to resample while applying an FIR filter to remove aliasing. Generally more accurate than simple resampling but slower than decimation.

Parameters:
  • new_sample_rate (float) – Target sample rate in samples per second.

  • pad_type (str, optional) – Padding method for edge effects: ‘mean’, ‘median’, ‘zero’ (default is ‘mean’).

  • inplace (bool, optional) – If True, modify current run. If False, return new run (default is False).

Returns:

If inplace is False, returns new resampled RunTS. Otherwise None.

Return type:

RunTS | None

Examples

Resample from 256 Hz to 100 Hz:

>>> resampled_run = run.resample_poly(100.0)
>>> print(resampled_run.sample_rate)
100.0

Notes

NaN values are filled with 0 before resampling. The polyphase method is particularly good for arbitrary sample rate ratios.

See also

decimate

Multi-stage decimation for downsampling

resample

Simple nearest-neighbor resampling

property run_metadata: Run

Run timeseries.

Returns:

Run-level metadata object (first run in first station).

Return type:

timeseries.Run

property sample_interval: float

Sample interval in seconds (inverse of sample_rate).

Returns:

Sample interval = 1 / sample_rate, or 0.0 if sample_rate is 0.

Return type:

float

Examples

>>> print(run.sample_interval)
0.00390625  # for 256 Hz
property sample_rate: float

Sample rate in samples per second.

Returns:

Sample rate estimated from time differences if data exists, otherwise from timeseries.

Return type:

float

Examples

>>> print(run.sample_rate)
256.0
set_dataset(array_list: list[ChannelTS] | list[DataArray] | Dataset, align_type: str = 'outer') None[source]

Set the dataset from a list of channels or existing dataset.

Creates an xarray Dataset from the input channels or dataset, validates metadata consistency, and updates dataset attributes with run metadata.

Parameters:
  • array_list (list[ChannelTS] | list[xr.DataArray] | xr.Dataset) – Input data as a list of ChannelTS objects, list of xarray DataArrays, or an existing xarray Dataset.

  • align_type (str, optional) –

    Method for aligning channels with different time indexes:

    • ’outer’ - use union of all time indexes (default)

    • ’inner’ - use intersection of time indexes

    • ’left’ - use time index from first channel

    • ’right’ - use time index from last channel

    • ’exact’ - raise ValueError if indexes don’t match exactly

    • ’override’ - rewrite indexes to match first channel (requires same size)

Notes

This method performs the following operations:

  1. Validates the input array_list

  2. Converts ChannelTS objects to xarray format

  3. Combines channels into a single Dataset

  4. Validates metadata consistency

  5. Updates dataset attributes with run metadata

When providing ChannelTS objects, their metadata and filters are automatically extracted and merged into the run’s metadata structure.

Examples

Set dataset from ChannelTS objects:

>>> ex = ChannelTS('electric', data=ex_data,
...                channel_metadata={'component': 'ex'})
>>> ey = ChannelTS('electric', data=ey_data,
...                channel_metadata={'component': 'ey'})
>>> run.set_dataset([ex, ey])

Set dataset with custom alignment:

>>> run.set_dataset([ex, ey, hx], align_type='inner')

Set dataset from existing xarray Dataset:

>>> import xarray as xr
>>> ds = xr.Dataset({'ex': ex_da, 'ey': ey_da})
>>> run.set_dataset(ds)

See also

dataset

Property for setting dataset with default alignment

add_channel

Add a single channel to existing dataset

_validate_array_list

Validation and conversion of array list

property start: MTime

Start time of the run in UTC.

Returns:

Start time from the dataset if data exists, otherwise from run_metadata.

Return type:

MTime

Examples

>>> print(run.start)
2020-01-01T00:00:00+00:00
property station_metadata: Station

Station timeseries.

Returns:

Station-level metadata object (first station in survey).

Return type:

timeseries.Station

property summarize_metadata: dict[str, any]

Get a summary of all channel timeseries.

Flattens the metadata from all channels into a single dictionary with keys in the format ‘channel.attribute’.

Returns:

Dictionary with flattened metadata from all channels.

Return type:

dict[str, any]

Examples

>>> meta_summary = run.summarize_metadata
>>> print(meta_summary.keys())
dict_keys(['ex.time_period.start', 'ex.sample_rate', ...])
property survey_metadata: Survey

Survey timeseries.

Returns:

Survey-level metadata object.

Return type:

timeseries.Survey

to_obspy_stream(network_code: str | None = None, encoding: str | None = None) Stream[source]

Convert time series to an ObsPy Stream object.

Creates an ObsPy Stream containing a Trace for each channel in the run.

Parameters:
  • network_code (str | None, optional) – Two-letter network code provided by FDSN DMC. If None, uses station timeseries.

  • encoding (str | None, optional) – Data encoding format (e.g., ‘STEIM2’, ‘FLOAT64’). If None, uses default encoding.

Returns:

Stream object containing Trace objects for all channels.

Return type:

obspy.core.Stream

Examples

>>> stream = run.to_obspy_stream(network_code='MT')
>>> print(stream)
3 Trace(s) in Stream:
MT.MT001..EX | 2020-01-01T00:00:00 - ... | 256.0 Hz, 4096 samples

See also

from_obspy_stream

Create RunTS from ObsPy Stream

ChannelTS.to_obspy_trace

Convert single channel

validate_metadata() None[source]

Check to make sure that the metadata matches what is in the data set.

updates metadata from the data.

Check the start and end times, channels recorded

wrangle_leap_seconds_from_obspy(array_list: list[ChannelTS]) list[ChannelTS][source]

Handle potential leap second issues from ObsPy streams.

Removes runs with only one sample that are numerically identical to adjacent samples, which may be artifacts of leap second handling.

Parameters:

array_list (list[ChannelTS]) – List of ChannelTS objects from ObsPy conversion.

Returns:

Filtered list with single-sample runs removed.

Return type:

list[ChannelTS]

Notes

This is experimental handling for issue #169. The exact behavior of ObsPy’s leap second handling is not fully documented.