Python signal.lfilter_zi方法代碼示例

# 需要導入模塊: from scipy import signal [as 別名]
# 或者: from scipy.signal import lfilter_zi [as 別名]
def compressor(x, thresh=-24, ratio=2, attackrel=0.045, sr=44100.0, dtype=np.float32):
    """
    simple compressor effect, code thanks to Eric Tarr @hackaudio
    Inputs:
       x:        the input waveform
       thresh:   threshold in dB
       ratio:    compression ratio
       attackrel:   attack & release time in seconds
       sr:       sample rate
    """
    attack = attackrel * sr  # convert to samples
    fc = 1.0/float(attack)     # this is like 1/attack time
    b, a = scipy_signal.butter(1, fc, analog=False, output='ba')
    zi = scipy_signal.lfilter_zi(b, a)

    dB = 20. * np.log10(np.abs(x) + 1e-6)
    in_env, _ = scipy_signal.lfilter(b, a, dB, zi=zi*dB[0])  # input envelope calculation
    out_env = np.copy(in_env)              # output envelope
    i = np.where(in_env >  thresh)          # compress where input env exceeds thresh
    out_env[i] = thresh + (in_env[i]-thresh)/ratio
    gain = np.power(10.0,(out_env-in_env)/20)
    y = x * gain
    return y 

# 需要導入模塊: from scipy import signal [as 別名]
# 或者: from scipy.signal import lfilter_zi [as 別名]
def geterrors(self, params):
        #copied from sandbox.tsa.arima.ARIMA
        p, q = self.nar, self.nma
        ar = np.concatenate(([1], -params[:p]))
        ma = np.concatenate(([1], params[p:p+q]))

        #lfilter_zi requires same length for ar and ma
        maxlag = 1+max(p,q)
        armax = np.zeros(maxlag)
        armax[:p+1] = ar
        mamax = np.zeros(maxlag)
        mamax[:q+1] = ma
        #remove zi again to match better with Skipper's version
        #zi = signal.lfilter_zi(armax, mamax)
        #errorsest = signal.lfilter(rhoy, rhoe, self.endog, zi=zi)[0] #zi is also returned
        errorsest = signal.lfilter(ar, ma, self.endog)
        return errorsest 

# 需要導入模塊: from scipy import signal [as 別名]
# 或者: from scipy.signal import lfilter_zi [as 別名]
def find_peak_vextors(price, return_ref=False, offest=0):
    """
    采用巴特沃斯信號濾波器，自適應尋找最佳極值點，決定平均周期分段數量。
    使用 scipy.Gaussian 機器學習統計算法進行第二次分析
    If you meet a Warning message, To slove this need upgrade scipy=>1.2. 
    but QUANTAXIS incompatible scipy=>1.2

    Parameters
    ----------
    price : (N,) array_like
        傳入需要查找極值點的價格-時間序列。
        The numerator coefficient vector of the filter.
    return_ref : bool or None, optional
        返回作為參照的平滑曲線，平滑曲線的目的是減少鋸齒抖動，減少計算的極值點。
        Return the smoothed line for reference.
    offest : int or None, optional
        傳遞參數時可能會被 .dropna() 或者 price[t:0] 等切片手段移除 nparray 頭部
        的 np.nan 元素，因為此函數返回的是向量節點的數組索引，為了跟原始參數對應，調用者
        可以指定一個補償偏移量，在返回的最大最小值中每個索引都會追加這個偏移量。
        The number of elements index offest, for jump np.nan in price's head.

    Returns
    -------
    x_tp_min, x_tp_max : ndarray
        包含最大值和最少值索引的數組
        The min/max peakpoint's index in array.

    """
    xn = price

    # Create an order 3 lowpass butterworth filter.
    b, a = butter(3, 0.05)

    # Apply the filter to xn.  Use lfilter_zi to choose the initial condition
    # of the filter.
    zi = lfilter_zi(b, a)
    z, _ = lfilter(b, a, xn, zi=zi * xn[0])

    # Apply the filter again, to have a result filtered at an order
    # the same as filtfilt.
    z2, _ = lfilter(b, a, z, zi=zi * z[0])

    # Use filtfilt to apply the filter.  If you meet a Warning need upgrade to
    # scipy=>1.2 but QUANTAXIS incompatible scipy=>1.2
    y = filtfilt(b, a, xn)

    # pass 1
    x_tp_min, x_tp_max = signal.argrelextrema(y, np.less)[0], signal.argrelextrema(y, np.greater)[0]
    n = int(len(price) / (len(x_tp_min) + len(x_tp_max))) * 2

    # peakutils 似乎一根筋隻能查最大極值，通過曲線反相的方式查找極小點
    mirrors = (price * -1) + np.mean(price) * 2

    # pass 2 使用 peakutils 查找
    x_tp_max = peakutils.indexes(price, thres=0.01 / max(price), min_dist=n)
    x_tp_min = peakutils.indexes(mirrors, thres=0.01 / max(price), min_dist=n)

    if (return_ref):
        return x_tp_min + offest, x_tp_max + offest, y
    else:
        return x_tp_min + offest, x_tp_max + offest