Skip to content

Adding code to read data and basic spectrum operations

Artem Basalaev requested to merge read_data_basic_spectrum_ops into main

Added core classes and organized basic workflow. It can be summarized as following:

  • with open function user specifies path to data files
  • DataFinder class takes care of finding files within and finding out data type (e.g. seismic time series data, or transmissibility spectrum) with the help of DataTypeClassifier
  • Empty SystemElement is created, the core class containing data and everything associated to an element in linear system. An element can represent, for example, a seismic noise input, or a passive isolation filter
  • Depending on data type, a data reader is used to read the data, currently either SeismicDataReader or TransmissibilityDataReader (both inherit from base class DataReader

At any time a spectrum can be generated from time series data using SystemElement.spectrum() which in turn either uses lpsd or returns the spectrum already calculated or loaded from external source.

A system element representing a signal/noise input can be multiplied by another one representing transmissibility/transfer function. This is achieved as follows:

  • Two system elements are multiplied with a * operator
  • Multiplication operator is overloaded in SystemElement. First it checks that there's one SystemElement with transfer function and another one with input signal/noise (order agnostic), otherwise multiplication cannot be performed
  • Multiplication itself is implemented in LinearMethods.multiply
  • For the transfer function spectrum, an interpolation with cubic splines is performed using a method in Spectrum class. It returns interpolated data points at a subset of frequencies of signal/noise that are within transfer function spectral range

Here is a use-case example:

import gravilab as gl
X1 = gl.open("seismic", channels = {"vertical":"_1_1.msd","north":"_1_2.msd","east":"_1_3.msd"})
X1.spectrum().data()["north"]["psd"].plot(logx=True, logy=True)
T = gl.open("transmissibility", channels = {"vertical":"1.","north":"2.","east":"3."})
T.spectrum().data()["north"].plot(logx=True, logy=True)
X2=X1*T
X2.spectrum().data()["north"]["psd"].plot(logx=True, logy=True)
Edited by Artem Basalaev

Merge request reports

Loading