Insert/append threshold crossing points (time and value) into a timeseries (pd.Series)

    Arguments:
      ts (pandas.Series): Time series of values to be interpolated at `thresh` crossings
      thresh (float or np.float64):

def insert_crossings(ts, thresh):
    """Insert/append threshold crossing points (time and value) into a timeseries (pd.Series)

    Arguments:
      ts (pandas.Series): Time series of values to be interpolated at `thresh` crossings
      thresh (float or np.float64):
    """
    # import time
    # tic0 = time.clock(); tic = tic0

    # int64 for fast processing, pandas.DatetimeIndex is 5-10x slower, 0.3 ms
    index = ts.index
    index_type = type(index)
    ts.index = ts.index.astype(np.int64)
    # toc = time.clock();
    # print((toc-tic)*1000); tic = time.clock()

    # value immediately before an upward thresh crossing, 6 ms
    preup = ts[(ts < thresh) & (ts.shift(-1) > thresh)]
    # toc = time.clock();
    # print((toc-tic)*1000); tic = time.clock()

    # values immediately after an upward thresh crossing, 4 ms\
    postup = ts[(ts.shift(1) < thresh) & (ts > thresh)]
    # toc = time.clock();
    # print((toc-tic)*1000); tic = time.clock()

    # value immediately after a downward thresh crossing, 1.8 ms
    postdown = ts[(ts < thresh) & (ts.shift(1) > thresh)]
    # toc = time.clock();
    # print((toc-tic)*1000); tic = time.clock()

    # value immediately before an upward thresh crossing, 1.9 ms
    predown = ts[(ts.shift(-1) < thresh) & (ts > thresh)]
    # toc = time.clock();
    # print((toc-tic)*1000); tic = time.clock()

    # upward slope (always positive) between preup and postup in units of
    # "value" per nanosecond (timestamps convert to floats as nanoseconds), 0.04 ms
    slopeup = (postup.values - preup.values) / (postup.index.values - preup.index.values).astype(np.float64)
    # toc = time.clock();
    # print((toc-tic)*1000); tic = time.clock()

    # upward crossing point index/time, 0.04 ms
    tup = preup.index.values + ((thresh - preup.values) / slopeup).astype(np.int64)
    # toc = time.clock();
    # print((toc-tic)*1000); tic = time.clock()

    # downward slope (always negative) between predown and postdown in units of
    # "value" per nanosecond (timestamps convert to floats as nanoseconds), 0.03 ms
    slopedown = (postdown.values - predown.values) / \
                (postdown.index.values - predown.index.values).astype(np.float64)
    # toc = time.clock();
    # print((toc-tic)*1000); tic = time.clock()

    # upward crossing point index/time, 0.02 ms
    tdown = predown.index.values + ((thresh - predown.values) / slopedown).astype(np.int64)
    # toc = time.clock();
    # print((toc-tic)*1000); tic = time.clock()

    # insert crossing points into time-series (if it had a regular sample period before, it won't now!), 2.0 ms
    ts.index = index  # pd.DatetimeIndex(ts.index)
    # toc = time.clock();
    # print((toc-tic)*1000); tic = time.clock()

    # insert crossing points into time-series (if it had a regular sample period before, it won't now!), 2.0 ms
    ts = ts.append(pd.Series(thresh * np.ones(len(tup)), index=index_type(tup.astype(np.int64))))
    # toc = time.clock();
    # print((toc-tic)*1000); tic = time.clock()

    # insert crossing points into time-series (if it had a regular sample period before, it won't now!), 1.9 ms
    ts = ts.append(pd.Series(thresh * np.ones(len(tdown)), index=index_type(tdown.astype(np.int64))))
    # toc = time.clock();
    # print((toc-tic)*1000); tic = time.clock()

    # if you don't `sort_index()`, numerical integrators in `scipy.integrate` will give the wrong answer, 0.1 ms
    ts = ts.sort_index()
    # toc = time.clock();
    # if you don't `sort_index()`, numerical integrators in `scipy.integrate` will give the wrong answer
    # print((toc-tic)*1000); tic = time.clock()
    # print((toc-tic0)*1000);
    return ts