удалить частоту питона (режекторный фильтр)? - PullRequest
Новая
       34

удалить частоту питона (режекторный фильтр)?

0 голосов
/ 09 мая 2018 
#complie by python3 only_test.py

import pyaudio
import numpy as np
import wave
import time
import math
#from pydub import AudioSegment
#from pydub.playback import play
#from scipy.signal import iirfilter
from scipy import signal

RATE = 48000
CHUNK = 4096
WIDTH = 2
volume = 0.0
duration = 1.0
#SHORT_NORMALIZE = (1.0/32768.0)
#INPUT_BLOCK_TIME = 1
#INPUT_BLOCK_PER_BLOCK = int(RATE*INPUT_BLOCK_TIME)

while True:

    #use a blackman window
    window = np.blackman(CHUNK)

    #load audio stream
    p = pyaudio.PyAudio()

    player = p.open(format=pyaudio.paInt16,
                    channels=1,
                    rate=RATE,
                    output=True,
                    frames_per_buffer=CHUNK)

    stream = p.open(format=pyaudio.paInt16,
                    channels=1,
                    rate=RATE,
                    input=True,
                    frames_per_buffer=CHUNK)

    #errorcount = 0

    for i in range(int(20*RATE/CHUNK)):

        sound = stream.read(CHUNK)

        #imp_ff = signal.filtfilt(b,a,sound)

        #playback microphone sound 
        #player.write(np.fromstring(sound,dtype=np.int16),CHUNK)

        #generate samples with return frequency to array
        #samples= (np.sin(2*np.pi*np.arange(RATE*duration)*freq/RATE)).astype(np.int16)

        #inverse frequency samples
        #inverse_samples = -samples

        #return frequency sound stream      
        #player.write(np.fromstring((volume*inverse_samples)\
        ,dtype=np.int16),CHUNK)

        #unpack the data and times by hamming window
        indata = np.array(wave.struct.unpack("%dh"%(len(sound)/WIDTH),\
                                             sound))*window

        #take the fft and square each value
        fftData = abs(np.fft.rfft(indata))*2

        #ifftData = abs(np.fft.irfft(indata))*2

        #find the maxium
        which = fftData[1:].argmax() + 1

        #use quadratic interpolation around the max
        if which != len(fftData)-1:
            y0,y1,y2 = np.log(fftData[which-1:which+2:])
            x1 = (y2-y0)*.5 / (2*y1-y2-y0)

            #find the frequency and output it
            freq = (which+x1)*RATE/CHUNK
            print("the freq is %d hz." % (freq))

        else:
            freq = which*RATE/CHUNK
            print("the freq is %d hz." % (freq))

        #playback the mic sound
        player.write(np.fromstring(sound,dtype=np.int16),CHUNK)

        if freq < 65:
           freq = 0

        #generate samples, note conversion to array
        #samples = 
        (np.sin(2*np.pi*np.arange(RATE*duration)*freq/RATE)).astype(np.int16)

        #invert phase of samples
        #result_samples = samples 

        #playback the invert_mic sound
        #player.write(np.fromstring(result_samples,dtype=np.int16),CHUNK)

    stream.stop_stream()
    stream.close()
    p.terminate()

В настоящее время мы обрабатываем микрофоны в режиме реального времени. Он предназначен для получения частоты через него и удаления синусоидального звука через режекторный фильтр (полосовой фильтр) для выходной частоты. Я не знаю, какой код написать, чтобы сделать режекторный фильтр (полосовой фильтр). У вас есть какой-нибудь код или библиотеки, чтобы помочь?

1 Ответ

0 голосов
/ 09 мая 2018

Поскольку вы уже используете scipy.signal, вы можете использовать scipy.signal.iirnotch . Возможно, вы также хотите прочитать некоторые сведения о БИХ-фильтрах и, например, Коэффициент качества .

Используйте его следующим образом: 

b, a = signal.iirnotch( w0, Q )

w0 - нормализованная частота.

Q - это коэффициент качества, который характеризует полосу пропускания режекторного фильтра -3 дБ относительно его центральной частоты.

Функция возвращает числитель b и знаменатель a полиномов БИХ-фильтра.

Пример: 

fs = 200.0  # Sample frequency (Hz)
f0 = 60.0  # Frequency to be removed from signal (Hz)
Q = 30.0  # Quality factor
w0 = f0 / (fs / 2 )  # Normalized Frequency
b, a = signal.iirnotch( w0, Q )
# Look at frequency response
w, h = signal.freqz( b, a )
freq = w * fs / ( 2 * np.pi )
plt.plot( freq, 20*np.log10( abs( h ) ) )
...