Source code for obspy.signal.differentiate_and_integrate

# -*- coding: utf-8 -*-
"""
Integration and differentiation routines.

:copyright:
    The ObsPy Development Team (devs@obspy.org)
:license:
    GNU Lesser General Public License, Version 3
    (https://www.gnu.org/copyleft/lesser.html)
"""
import numpy as np
import scipy.integrate
import scipy.interpolate


[docs] def integrate_cumtrapz(data, dx, **kwargs): """ Performs first order integration of data using the trapezoidal rule. :param data: Data array to integrate. :param dx: Sample spacing usually in seconds. """ # Integrate. Set first value to zero to avoid changing the total # length of the array. # (manually adding the zero and not using `cumtrapz(..., initial=0)` is a # backwards compatibility fix for scipy versions < 0.11. ret = scipy.integrate.cumulative_trapezoid(data, dx=dx) return np.concatenate([np.array([0], dtype=ret.dtype), ret])
[docs] def integrate_spline(data, dx, k=3, **kwargs): """ Integrate by generating an interpolating spline and integrating that. :param data: The data to integrate. :param dx: Sample spacing usually in seconds. :param k: Spline order. 1 is linear, 2 quadratic, 3 cubic, Must be between 1 and 5. """ time_array = np.linspace(0, (len(data) - 1) * dx, len(data)) spline = scipy.interpolate.InterpolatedUnivariateSpline(time_array, data, k=k) # Backport of the antiderivative() method for scipy versions < 0.13.0. # Can be removed once the minimum supported version is equal or larger # to this. if not hasattr(spline, "antiderivative"): t, c, k = spline._eval_args # Compute the multiplier in the antiderivative formula. dt = t[k + 1:] - t[:-k - 1] # Compute the new coefficients c = np.cumsum(c[:-k - 1] * dt) / (k + 1) c = np.r_[0, c, [c[-1]] * (k + 2)] # New knots t = np.r_[t[0], t, t[-1]] k += 1 tmp = scipy.interpolate.InterpolatedUnivariateSpline.__new__( scipy.interpolate.InterpolatedUnivariateSpline) tmp._eval_args = t, c, k tmp._data = (None, None, None, None, None, k, None, len(t), t, c, None, None, None, None) tmp.ext = 0 return tmp(time_array) return spline.antiderivative(n=1)(time_array)