processing.py 文件源码

def get_multievent_sg(cum_trace, tmin, tmax, tstart,
                      threshold_inside_tmin_tmax_percent,
                      threshold_inside_tmin_tmax_sec, threshold_after_tmax_percent):
    """
        Returns the tuple (or a list of tuples, if the first argument is a stream) of the
        values (score, UTCDateTime of arrival)
        where scores is: 0: no double event, 1: double event inside tmin_tmax,
            2: double event after tmax, 3: both double event previously defined are detected
        If score is 2 or 3, the second argument is the UTCDateTime denoting the occurrence of the
        first sample triggering the double event after tmax
        :param trace: the input obspy.core.Trace
    """
    tmin = utcdatetime(tmin)
    tmax = utcdatetime(tmax)
    tstart = utcdatetime(tstart)

    # split traces between tmin and tmax and after tmax
    traces = [cum_trace.slice(tmin, tmax), cum_trace.slice(tmax, None)]

    # calculate second derivative and normalize:
    second_derivs = []
    max_ = np.nan
    for ttt in traces:
        ttt.taper(type='cosine', max_percentage=0.05)
        sec_der = savitzky_golay(ttt.data, 31, 2, deriv=2)
        sec_der_abs = np.abs(sec_der)
        idx = np.nanargmax(sec_der_abs)
        # get max (global) for normalization:
        max_ = np.nanmax([max_, sec_der_abs[idx]])
        second_derivs.append(sec_der_abs)

    # normalize second derivatives:
    for der in second_derivs:
        der /= max_

    result = 0

    # case A: see if after tmax we exceed a threshold
    indices = np.where(second_derivs[1] >= threshold_after_tmax_percent)[0]
    if len(indices):
        result = 2

    # case B: see if inside tmin tmax we exceed a threshold, and in case check the duration
    deltatime = 0
    indices = np.where(second_derivs[0] >= threshold_inside_tmin_tmax_percent)[0]
    starttime = endtime = None
    if len(indices) >= 2:
        idx0 = indices[0]
        idx1 = indices[-1]
        starttime = timeof(traces[0], idx0)
        endtime = timeof(traces[0], idx1)
        deltatime = endtime - starttime
        if deltatime >= threshold_inside_tmin_tmax_sec:
            result += 1

    return result, deltatime, starttime, endtime