Python signal.ellip函数代码示例

    def __amp_detect(self, x):
        
        ref = np.floor(self.min_ref_per*self.sr/1000.0)
        
        # HIGH-PASS FILTER OF THE DATA
        (b,a) = signal.ellip(2, 0.1, 40, [self.fmin_detect*2.0/self.sr,self.fmax_detect*2.0/self.sr], btype='bandpass', analog=0, output='ba')
        xf_detect = signal.filtfilt(b, a, x)
        (b,a) = signal.ellip(2, 0.1, 40, [self.fmin_sort*2.0/self.sr,self.fmax_sort*2.0/self.sr], btype='bandpass', analog=0, output='ba')
        xf = signal.filtfilt(b, a, x)
        
        
        noise_std_detect = scipy.median(np.abs(xf_detect))/0.6745;
        noise_std_sorted = scipy.median(np.abs(xf))/0.6745;
       
        thr = self.stdmin * noise_std_detect        #thr for detection is based on detected settings.
        thrmax = self.stdmax * noise_std_sorted     #thrmax for artifact removal is based on sorted settings.
        
        # LOCATE SPIKE TIMES
        nspk = 0;
        xaux = np.argwhere(xf_detect[self.w_pre+1:len(xf_detect)-self.w_post-1-1] > thr) + self.w_pre + 1
        xaux = np.resize(xaux,len(xaux))
        xaux0 = 0;
        index = []
        for i in range(len(xaux)):
            if xaux[i] >= (xaux0 + ref):
            # after find a peak it begin search after ref over the last xaux
                iaux = xf[xaux[i]:xaux[i]+np.floor(ref/2.0)].argmax(0)    # introduces alignment
                nspk = nspk + 1
                index.append(iaux + xaux[i])
                xaux0 = index[nspk-1];
        
        # SPIKE STORING (with or without interpolation)
        ls = self.w_pre + self.w_post
        spikes = np.zeros([nspk,ls+4])
        xf = np.concatenate((xf,np.zeros(self.w_post)),axis=0)
        
        for i in range(nspk):                          # Eliminates artifacts
            if np.max( np.abs( xf[index[i]-self.w_pre:index[i]+self.w_post] )) < thrmax :
                spikes[i,:] = xf[index[i]-self.w_pre-1:index[i]+self.w_post+3]
     
        aux = np.argwhere(spikes[:,self.w_pre] == 0)       #erases indexes that were artifacts
        if len(aux) != 0:
            aux = aux.reshape((1,len(aux)))[0]
            spikes = np.delete(spikes, aux, axis = 0)
            index = np.delete(index,aux)
 
        if self.interpolation == 'y':
            # Does interpolation
            spikes = self.__int_spikes(spikes)

        return spikes, thr, index

 def __init__(self, timestep):
     self.sampling_rate = int(1. / timestep)
     self.timestep = timestep
     self.c_detect = ellip(2, .1, 40,
                           (2 * timestep * DETECT_LOW,
                            2 * timestep * DETECT_HIGH),
                           'bandpass')
     self.c_extract = ellip(2, .1, 40,
                            (2 * timestep * EXTRACT_LOW,
                             2 * timestep * EXTRACT_HIGH),
                            'bandpass')
     self.c_notch = ellip(2, .5, 20,
                          (2 * timestep * 1999, 2 * timestep * 2001),
                          'bandstop')

 def HPman(self, fil_dict):
     """Elliptic HP filter, manual order"""
     self._get_params(fil_dict)
     if not self._test_N():
         return -1
     self._save(fil_dict, sig.ellip(self.N, self.A_PB, self.A_SB, self.F_PB,
                     btype='highpass', analog=self.analog, output=self.FRMT))

    def _sos(self, sfreq):
        nyq = sfreq / 2.
        low_stop, low_pass, high_pass, high_stop, gpass, gstop = self.args
        if high_stop is None:
            assert low_stop is not None
            assert high_pass is None
        else:
            high_stop /= nyq
            high_pass /= nyq

        if low_stop is None:
            assert low_pass is None
        else:
            low_pass /= nyq
            low_stop /= nyq

        if low_stop is None:
            btype = 'lowpass'
            wp, ws = high_pass, high_stop
        elif high_stop is None:
            btype = 'highpass'
            wp, ws = low_pass, low_stop
        else:
            btype = 'bandpass'
            wp, ws = (low_pass, high_pass), (low_stop, high_stop)
        order, wn = signal.ellipord(wp, ws, gpass, gstop)
        return signal.ellip(order, gpass, gstop, wn, btype, output='sos')

	def bp_filt(self, df):
		'''bandpass filter'''
		cutfreq = np.array([self.up_Fcut, self.low_Fcut])
		b, a = signal.ellip(4, 0.1, 40, 
							cutfreq*2/self.dataset.Fs, btype='bandpass')
		output = self.filtfilt_df(df, a, b)
		return output

 def LPmin(self, fil_dict):
     """Elliptic LP filter, minimum order"""
     self.get_params(fil_dict)
     self.N, self.F_PBC = ellipord(self.F_PB,self.F_SB, self.A_PB,self.A_SB,
                                                       analog = self.analog)
     self.save(fil_dict, sig.ellip(self.N, self.A_PB, self.A_SB, self.F_PBC,
                         btype='low', analog = self.analog, output = frmt))

 def design_filter(self, fs):
     if np.isinf(self.LPFcutoff):
         N, Ws = ellipord(self.HPFcutoff * 1e3 / fs * 2,
                          max(5 * 1e3 / fs * 2,
                              (self.HPFcutoff - 5) * 1e3 / fs * 2),
                          self.Rp,
                          self.Rs)
         b, a = ellip(N, self.Rp, self.Rs, Ws, 'high')
     else:
         N, Ws = ellipord([self.HPFcutoff * 1e3 / fs * 2, self.LPFcutoff * 1e3 / fs * 2],
                          [max(5*1e3/fs*2,(self.HPFcutoff-5)*1e3/fs*2),
                           min((fs/2-5e3)/fs*2,(self.LPFcutoff+5)*1e3/fs*2)],
                          self.Rp,
                          self.Rs)
         b, a = ellip(N, self.Rp, self.Rs, Ws)
     return b, a

def pitch_filter_bank(ratios_pitches=None, fs=16000, Q=25.0, max_loss_pass=1.0, min_attenuation_stop=50.0):
    """
    lowest pitch: 20.6 Hz = pitch 16, the lowest pitch above the low threshold of hearing
    highest pitch: 7458.6 Hz = pitch 118, the highest pitch below half of the sampling frequency (fs = 16000Hz)
    Note that 119 is technically below the nyquist frequency (~7900Hz), but the right stopband frequency wouldn't be.
    
    fs: sampling frequency of the input in Hz
    Q: Q factor = frequency / bandwidth, used to determine passband and stopband frequencies of the elliptic filters
    max_loss_pass: maximal loss in passband in dB
    min_attenuation_stop: minimal attenuation in stopband in dB
    """
    if ratios_pitches is None:
        ratios_pitches = RATIOS_PITCHES_DEFAULT
        # structure: tuples of sampling frequency ratios and sets of pitches

    filters = {} # dictionary indexed by sampling frequency ratio. Each item is again a dictionary indexed by pitch, giving a filter coefficient tuple.
    
    for ratio, pitches in ratios_pitches:
        filters[ratio] = {}
        current_fs = float(fs / ratio) # input sampling frequency for the current set of pitches
        nyquist_freq = current_fs / 2
        for pitch in pitches:
            freq = pitch2freq(pitch)
            w = freq / nyquist_freq # omega = normalised frequency
            w_pass = (w * (1 - 1 / (2*Q)), w * (1 + 1 / (2*Q)))
            w_stop = (w * (1 - 1 / Q), w * (1 + 1 / Q))
            n, w_natural = sig.ellipord(w_pass, w_stop, max_loss_pass, min_attenuation_stop)
            coeff_b, coeff_a = sig.ellip(n, max_loss_pass, min_attenuation_stop, w_natural, btype='bandpass') # get filter coefficients
            # note that scipy's ellip differs from matlab's in that it will always generate a lowpass filter by default.
            # btype='bandpass' needs to be passed explicitly!
            filters[ratio][pitch] = (coeff_b, coeff_a)
    
    return filters

 def BPmin(self, fil_dict):
     """Elliptic BP filter, minimum order"""
     self.get_params(fil_dict)
     self.N, self.F_PBC = ellipord([self.F_PB, self.F_PB2],
         [self.F_SB, self.F_SB2], self.A_PB, self.A_SB, analog = self.analog)
     self.save(fil_dict, sig.ellip(self.N, self.A_PB, self.A_SB, self.F_PBC,
                     btype='bandpass', analog = self.analog, output = frmt))

    def data_hpass(self, x, Wp, srate):
        ''' High-pass filter '''
        Wp = float(Wp*2/srate)
        Ws = Wp*float(self.lineEdit_19.text())
        Rp = float(self.lineEdit_17.text())
        Rs = float(self.lineEdit_18.text())

        tempstring = self.lineEdit_16.text()
        if tempstring == 'auto':
            if self.comboBox_2.currentIndex() == 0:
                (norder, Wn) = buttord(Wp, Ws, Rp, Rs)
            elif self.comboBox_2.currentIndex() == 1:
                (norder, Wn) = ellipord(Wp, Ws, Rp, Rs)
            else:
                (norder, Wn) = cheb1ord(Wp, Ws, Rp, Rs)
        else:
            norder = float(tempstring)
            Wn = Wp

        if self.comboBox_2.currentIndex() == 0:
            (b, a)  =  butter(norder, Wn, btype = 'high')
        elif self.comboBox_2.currentIndex() == 1:
            (b, a)  =  ellip(norder, Rp, Rs, Wn)
        else:
            (b, a)  =  cheby1(norder, Rp, Wn)


        y  =  filtfilt(b, a, x)

        return(y)

    def __init__(self, type, fc, gain = 0, Q = 1, enabled = True):
        self._enabled = enabled
        self._type = type
        self._fc = fc
        self._g = gain
        self._Q = Q

        if type == FilterType.HPBrickwall:
            z, p, k = scsig.ellip(12, 0.01, 80, fc, 'high', output='zpk')
            self._sos = zpk2sos(z, p, k)[0]
        elif type == FilterType.LPBrickwall:
            z, p, k = scsig.ellip(12, 0.01, 80, fc, 'low', output='zpk')
            self._sos = zpk2sos(z, p, k)[0]
        elif type == FilterType.HPButter:
            z, p, k = scsig.butter(2 ** Q, fc, btype = 'high', output='zpk')
            self._sos = zpk2sos(z, p, k)[0]
        elif type == FilterType.LPButter:
            z, p, k = scsig.butter(2 ** Q, fc, output='zpk')
            self._sos = zpk2sos(z, p, k)[0]
        elif type == FilterType.LShelving or type == FilterType.HShelving:
            A = 10 ** (gain / 20)
            wc = np.pi * fc
            wS = np.sin(wc)
            wC = np.cos(wc)
            alpha = wS / (2 * Q)
            beta = A ** 0.5 / Q
            c = 1
            if type == FilterType.LShelving:
                c = -1

            b0 = A * (A + 1 + c * (A - 1) * wC + beta * wS)
            b1 = - c * 2 * A * (A - 1 + c * (A + 1) * wC)
            b2 = A * (A + 1 + c * (A - 1) * wC - beta * wS)
            a0 = (A + 1 - c * (A - 1) * wC + beta * wS)
            a1 = c * 2 * (A - 1 - c * (A + 1) * wC)
            a2 = (A + 1 - c * (A - 1) * wC - beta * wS)
            self._sos = np.array([[ b0, b1, b2, a0, a1, a2 ]])
        elif type == FilterType.Peak:
            self.g = gain
            wc = np.pi * fc
            b, a = scsig.bilinear([1, 10 ** (gain / 20) * wc / Q, wc ** 2],
                [1, wc / Q, wc ** 2])
            self._sos = np.append(b, a).reshape(1, 6)

        self._ord = self._sos.shape[0] * 2
        self.icReset()

 def BSman(self, fil_dict):
     """Elliptic BS filter, manual order"""
     self._get_params(fil_dict)
     if not self._test_N():
         return -1
     self._save(fil_dict, sig.ellip(self.N, self.A_PB, self.A_SB, 
         [self.F_PB,self.F_PB2], btype='bandstop', analog=self.analog, 
                                                             output=self.FRMT))

 def BPmin(self, fil_dict):
     """Elliptic BP filter, minimum order"""
     self._get_params(fil_dict)
     self.N, self.F_PBC = ellipord([self.F_PB, self.F_PB2],
         [self.F_SB, self.F_SB2], self.A_PB, self.A_SB, analog=self.analog)
     if not self._test_N():
         return -1
     self._save(fil_dict, sig.ellip(self.N, self.A_PB, self.A_SB, self.F_PBC,
                     btype='bandpass', analog=self.analog, output=self.FRMT))

    def _design(self):
        if not self.stopband_attenuation:
            self.stopband_attenuation = self.filter_parameters['stopband_attenuation']

        if not self.ripple:
            self.ripple = self.filter_parameters['ripple']

        if self.already_normalized_Wn:
            self.Z, self.P, self.K = signal.ellip(self.N, self.ripple,
                                                  self.stopband_attenuation,
                                                  self.Wn,
                                                  self.filter_kind,
                                                  analog=False, output='zpk')
        else:
            self.Z, self.P, self.K = signal.ellip(self.N, self.ripple,
                                                  self.stopband_attenuation,
                                                  self.normalize_Wn(),
                                                  self.filter_kind, analog=False,
                                                  output='zpk')

    def HPmin(self, fil_dict):
        """Elliptic HP filter, minimum order"""
        self._get_params(fil_dict)
        self.N, self.F_PBC = ellipord(self.F_PB,self.F_SB, self.A_PB,self.A_SB,
                                                          analog=self.analog)
#       force even N
        if (self.N%2)== 1:
            self.N += 1
        self._save(fil_dict, sig.ellip(self.N, self.A_PB, self.A_SB, self.F_PBC,
                        btype='highpass', analog=self.analog, output=self.FRMT))