obspy.realtime.rttrace.RtTrace

class RtTrace(max_length=None, *args, **kwargs)[source]

Bases: obspy.core.trace.Trace

An object containing data of a continuous series constructed dynamically from sequential data packets.

New data packets may be periodically appended. Registered time-domain processes can be applied to the new data and the resulting trace will be left trimmed to maintain a specified maximum trace length.

Parameters:max_length (int, optional) maximum trace length in seconds

Example

RtTrace has been built to handle real time processing of periodically append data packets, such as adding and processing data requested from an SeedLink server. See obspy.seedlink for further information.

For the sake of simplicity we will just split data of an existing example file into multiple chucks (Trace objects) of about equal size (step 1 + 2) and append those chunks in a simple loop (step 4) into an RtTrace object. Additionally there are two real time processing functions registered to the RtTrace object (step 3) which will automatically process any appended data chunks.

  1. Read first trace of example SAC data file and extract contained time offset and epicentral distance of an earthquake:

    >>> import numpy as np
    >>> from obspy.realtime import RtTrace
    >>> from obspy import read
    >>> from obspy.realtime.signal import calculateMwpMag
    >>> data_trace = read('/path/to/II.TLY.BHZ.SAC')[0]
    >>> len(data_trace)
    12684
    >>> ref_time_offset = data_trace.stats.sac.a
    >>> print(ref_time_offset)
    301.506
    >>> epicentral_distance = data_trace.stats.sac.gcarc
    >>> print(epicentral_distance)
    30.0855
    
  2. Split given trace into a list of three sub-traces:

    >>> traces = data_trace / 3
    >>> [len(tr) for tr in traces]
    [4228, 4228, 4228]
    
  3. Assemble real time trace and register two processes:

    >>> rt_trace = RtTrace()
    >>> rt_trace.registerRtProcess('integrate')
    1
    >>> rt_trace.registerRtProcess('mwpIntegral', mem_time=240,
    ...     ref_time=(data_trace.stats.starttime + ref_time_offset),
    ...     max_time=120, gain=1.610210e+09)
    2
    
  4. Append and auto-process packet data into RtTrace:

    >>> for tr in traces:
    ...     processed_trace = rt_trace.append(tr, gap_overlap_check=True)
    ...
    >>> len(rt_trace)
    12684
    
  5. Some post processing to get Mwp:

    >>> peak = np.amax(np.abs(rt_trace.data))
    >>> print(peak)
    0.136404
    >>> mwp = calculateMwpMag(peak, epicentral_distance)
    >>> print(mwp)
    8.78902911791
    

Attributes

__dict__
__doc__ str(object) -> string
__module__ str(object) -> string
__weakref__ list of weak references to the object (if defined)
have_appended_data bool(x) -> bool
id Returns a SEED compatible identifier of the trace.

Public Methods

append Appends a Trace object to this RtTrace.
copy Returns a deepcopy of this RtTrace.
count Returns number of data samples of the current trace.
decimate Downsample trace data by an integer factor.
detrend Method to remove a linear trend from the trace.
differentiate Method to differentiate the trace with respect to time.
filter Filters the data of the current trace.
getId Returns a SEED compatible identifier of the trace.
integrate Method to integrate the trace with respect to time.
max Returns the value of the absolute maximum amplitude in the trace.
normalize Method to normalize the trace to its absolute maximum.
plot Creates a simple graph of the current trace.
registerRtProcess Adds real-time processing algorithm to processing list of this RtTrace.
resample Resample trace data using Fourier method.
rtProcessFunctionsToString Return doc string for all predefined real-time processing functions.
simulate Correct for instrument response / Simulate new instrument response.
slice Returns a new Trace object with data going from start to end time.
spectrogram Creates a spectrogram plot of the trace.
split Splits Trace object containing gaps using a NumPy masked array into
std Method to get the standard deviation of amplitudes in the trace.
taper Method to taper the trace.
times For convenient plotting compute a Numpy array of seconds since
trigger Runs a triggering algorithm on the data of the current trace.
trim Cuts current trace to given start and end time.
verify Verifies current trace object against available meta data.
write Saves current trace into a file.

Private Methods

_addProcessingInfo Adds the given informational string to the processing field in the
_ltrim Cuts current trace to given start time.
_rtrim Cuts current trace to given end time.

Special Methods

__add__ Too ambiguous, throw an Error.
__div__ Splits Trace into new Stream containing num Traces of the same size.
__eq__ Implements rich comparison of RtTrace objects for “==” operator.
__ge__ Too ambiguous, throw an Error.
__getitem__ __getitem__ method of Trace object.
__gt__ Too ambiguous, throw an Error.
__init__ Initializes an RtTrace.
__le__ Too ambiguous, throw an Error.
__len__ Returns number of data samples of the current trace.
__lt__ Too ambiguous, throw an Error.
__mod__ Splits Trace into new Stream containing Traces with num samples.
__mul__ Creates a new Stream containing num copies of this trace.
__ne__ Implements rich comparison of Trace objects for ”!=” operator.
__original_str__ Returns short summary string of the current trace.
__setattr__ __setattr__ method of Trace object.
__str__ Monkey patch for the __str__ method of the Trace object. SEGY object do not

This Page