1

My question is related to the explanation here by A. Levy: Analyze audio using Fast Fourier Transform

How can I produce a bandpass filter on these complex numbers... 

[-636.00000000 +0.00000000e+00j  -47.84161618 -2.80509841e+02j
30.69754505 -1.30624718e+01j -109.94022791 +7.58155488e+00j
-3.18538186 +1.44880882e+01j -120.36687555 +5.45225425e+00j
50.48671763 +1.69504204e+01j   31.56751791 -7.22728042e+01j
-17.96079093 -3.17853727e+01j  -19.25527276 +5.08151876e+00j
18.38143611 -2.60879726e+01j  -27.15617871 -4.39789289e+01j...

... knowing already that I will then convert my complex array back into the time domain using ifft. My band pass should only allow 18500Hz to 19500Hz through (super-sonic, I know. I'm specifically looking for those tones to arrive from a tone generator across the room). My sampling rate it 44100Hz and the number of points in my FFT is 128.

I will then be using the frequency detection code here (which already works on my non-filterd samples) Python frequency detection in order to give me the pitch of the audio tone when it arrives. If there is no tone played at 19000HZ (for example), then the output of the pitch detection function should return 0.

My current filter code (which is not giving me the correct result):

samp = [-10 -16 -15 -11 -8 -10 -9  -12 -7 -13 -4 -10 
        -1  -4  -11 -6 -8 -8 -10 -6 -9 -7 -16 -11 5 -14   
        -9 -3 -9 -7 -7 -6 -3 -11 -13 -9 -10 -4 -6 -7 
        -11 -15 -15 -3 -5 -15 -11 -8 -13 -9 -12 -10 -8 -16  
        -13 -5 -4 -12 -14 -8 -14 -6 -7 -7 -4 -6 -9 -4 -4  
        -4  1 -10 -3 -9 -9  -1 -5 -2 -5 -3 -3  2 -3  2  
        -5 -4 -6  1 -2 -6 -8  -3 -10 -11 -6 -2 -5  -3   
         0  3  0  1  1 -1 -3 -3  1  3 -3 -3  3 -3 -1  
        -3  -1  2  1  0 -8  0  6 -3 -4 -7 -5 -10 -4 -4]

sample_time = 0.000022675   # 1/44100

low_cut = 18500
high_cut = 19500

float_samp = np.float32(samp)
fft_spectrum = np.fft.fft(float_samp)

freqs = np.fft.fftfreq(len(fft_spectrum), d=sample_time)

### wrong approach ###
for i in range(len(freqs)): # (H-red)

    if abs(freqs[i]) >= high or abs(freqs[i]) <= low:
        fft_spectrum[i] = 0.0 
### -------------- ###

time_domain = np.fft.ifft(fft_spectrum)
converted = np.int16(time_domain)

It's currently only returning values between 18500Hz and 19500Hz after pitch detection, whether or not the tone is present. I believe this is the case because I'm knocking out all the information from my fft_spectrum list except those within with my filter values. Even though the tone is not there, the only information left is the within that band, so that's what the pitch detector is reading. This is my assumption.

Please Note: I can not use Scipy as I'm deploying this on Android and that library is not available on that platform. I'm hoping there is a way to do this through Numpy.

Community
  • 1
  • 1
A. Wylie - DigiCoyote
  • 67
  • 1
  • 13

1 Answers1

1

The frequency detection code you linked is performing an FFT and then finding the bin with the largest magnitude. The only way it's going to return a zero is if the largest magnitude is in the 0th bin. Even though you may not be producing a tone at the frequency you are interested in there is most certainly energy there and therefore is a candidate for being the biggest. To do you what you're asking you'll need to modify the referenced code to provide the desired behavior. For example you could apply some kind of minimum threshold.

which = fftData[1:].argmax() + 1
if (fftData[which] < threshold)
    return -1; // no peak found

While you're at it: I don't understand why you are going through the process of bandpass filtering the signal when you could just limit the frequency detection to search the bins of interest:

binMin = floor(low_cut / 22050.0 * 128)
binMax = ceil(high_cut / 22050.0 * 128)
which = fftData[binMin:binMax].argmax() + 1
jaket
  • 9,140
  • 2
  • 25
  • 44
  • I like where you're headed with this. 1.) If I just take the bin approach, being your second piece of code, would I need to set the threshold at all like you showed me with your first set of code? I can see this app needing to run in a very loud environment, and being able to detect those high frequencies above the regular room volume (which may or may not be louder that the super-sonic tones themselves). – A. Wylie - DigiCoyote May 05 '15 at 03:47
  • 2.) Are binMin and binMax supposed to equal 0 and 1 respectively? That's what they seem to be coming out as, and none of the frequencies seem to be getting filtered. – A. Wylie - DigiCoyote May 05 '15 at 03:54
  • Yes you would still need to apply some kind of threshold since even though the numbers may still be very small if you simply ask for the max your going to get something. One approach would be to look at the bins with no signal to get an idea of the noise floor and then look for your signal to be sufficiently above that. – jaket May 05 '15 at 03:54
  • oh yeah there's a mistake there. I'll edit the answer. – jaket May 05 '15 at 03:57
  • Maybe I need to do: binMax = math.ceil(high_cut / 22050.0 * 128.0 / 2)... I think I need to divide by 2 because fftData is 65 and samp is 128, yes? Otherwise it's throwing a *ValueError: attempt to get argmax of an empty sequence* – A. Wylie - DigiCoyote May 05 '15 at 04:52
  • Ah, it's because I'm doing 'rfft' instead of 'fft'. I don't need to divide by two if I use 'fft'. The original code uses 'rfft' though. – A. Wylie - DigiCoyote May 05 '15 at 04:56
  • do you mean 64? I read earlier in your post that the # of points in your fft was 128. If you input 128 samples and get out 64 bins then you need to divide by 64. Hope that helps. – jaket May 05 '15 at 04:57
  • also, you should be applying a window to the time-domain signal. – jaket May 05 '15 at 04:59
  • It totally helps. I know it should be 64, but pycharm's debugger says that the ndarray is 65. Odd. Now that I've wrapped my head around this a bit, potentially I don't need to be doing any pitch detection, I'm just looking for changes in power of the frequency bins I'm interested in? – A. Wylie - DigiCoyote May 05 '15 at 05:04
  • Looks like my magnitude tests are coming back positive. Thanks for pointing me in that direction. :) Basically I'm testing each bin's magnitude against the one in the previous sample and when I find an increase, I know my signal has arrived. You had mentioned putting a window on the time-domain signal before I fft it. Without Scipy, how would you recommend going about doing that? (Id like to move this conversation to chat, but I need 20 rep. I'm new here) – A. Wylie - DigiCoyote May 05 '15 at 06:57
  • Let us [continue this discussion in chat](http://chat.stackoverflow.com/rooms/77059/discussion-between-a-wylie-digicoyote-and-jaket). – A. Wylie - DigiCoyote May 06 '15 at 05:51