Python特定频率删除(陷波滤波器)? [英] python specific frequency remove(notch filter)?

问题描述

#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()

我们目前正在实时处理麦克风. 它旨在通过该频率获得频率，并通过陷波滤波器(带阻滤波器)消除输出频率的正弦波声音. 我不知道写什么代码来做陷波滤波器(带阻滤波器). 您有任何代码或库可以帮助您吗?

We are currently processing microphones in real time. It is designed to obtain the frequency through it and to remove the sine wave sound through the notch filter (bandstop filter) for the output frequency. I do not know what code to write to do a notch filter (bandstop filter). Do you have any code or libraries to help?

推荐答案

虽然Arude的解决方案很不错，但我还是简化了一下.另外，在scipy docs中，w0术语可能会造成混乱.它说必须对它进行归一化，但是只要您输入采样率，它就会为您完成.

While ARude's solution is decent, I made this one to be a bit more simple. Also, in the scipy docs the w0 term can be pretty confusing. It says it must be normalized, but it does it for you if you just enter the sampling rate.

我将这种简化形式包装到以下代码中，并添加了一个带有虚假信号的完整示例，其中我消除了60Hz的嗡嗡声:

I've wrapped such simplifications into the following code, and added a full working example with a fake signal where I remove the 60Hz hum:

from scipy import signal
import matplotlib.pyplot as plt
import numpy as np

# Create/view notch filter
samp_freq = 1000  # Sample frequency (Hz)
notch_freq = 60.0  # Frequency to be removed from signal (Hz)
quality_factor = 30.0  # Quality factor
b_notch, a_notch = signal.iirnotch(notch_freq, quality_factor, samp_freq)
freq, h = signal.freqz(b_notch, a_notch, fs = samp_freq)
plt.figure('filter')
plt.plot( freq, 20*np.log10(abs(h)))

# Create/view signal that is a mixture of two frequencies
f1 = 17
f2 = 60
t = np.linspace(0.0, 1, 1_000)
y_pure = np.sin(f1 * 2.0*np.pi*t) + np.sin(f2 * 2.0*np.pi*t) 
plt.figure('result')
plt.subplot(211)
plt.plot(t, y_pure, color = 'r')

# apply notch filter to signal
y_notched = signal.filtfilt(b_notch, a_notch, y_pure)

# plot notch-filtered version of signal
plt.subplot(212)
plt.plot(t, y_notched, color = 'r')

这篇关于Python特定频率删除(陷波滤波器)?的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！