5

I'm really new on sound processing, so maybe my question will be trivial. What I want to do is to extract a specific frequency range (let's say 150-400 Hz) from a wav file, using R. In other words, I want to create another wave file (wave2) that contains only the frequency component that I specify (150 to 400 Hz, or what else).

I read something on the net, and I discovered out that this can be done with a FFT analysis, and here's come the problems.

Suppose I've this code:

library(sound)
s1 <- Sine(440, 1)
s2 <- Sine(880, 1)
s3 <- s1 + s2

s3.s <- as.vector(s3$sound)
  # s3.s is now a vector, with length 44100; 
  # bitrate is 44100 (by default)
  # so total time of s3 is 1sec.

  # now I calculate frequencies
N <- length(s3.s)   # 44100
k <- c(0:(N-1))
Fs <- 44100         # sampling rate
T <- N / Fs
freq <- k / T
x <- fft(s3.s) / N

plot(freq[1:22050], x[1:22050], type="l") # we need just the first half of FFT computation

The plot we obtain is:

enter image description here

Well, there are two peaks. If we want to know to what frequency they correspond, just find:

order(Mod(x)[1:22050], decreasing=T)[1:10]
[1] 441 881 882 880 883 442 440 879 884 878

First two values are really near to the frequency I've used to create my sound:

        real     computed
 Freq1: 440   |  441 
 Freq2: 880   |  881

So, now comes the problem: how to proceed, if I want to delete from my sound the frequencies in the range, say, (1, 500) ? And how to select (and save) only the range (1, 500) ? What I attend, is that my new sound (with deleted frequencies) will be something near to simple Sine(freq=880, duration=1) (I know, it cannot be exactly like so!). Is that possible?

I'm pretty sure that fft(DATA, inverse = TRUE) is what I need. But I'm not sure, and however I don't know how to proceed.

Tommaso
  • 527
  • 1
  • 5
  • 17

3 Answers3

3

If you don't want to programm it, you can use Praat.

Praat is a free scientific software program for the analysis of speech in phonetics. But you can also use it to edit the spectrum of any sound (remove frequencies, ...) and then export the result as a new sound file.

user1493046
  • 352
  • 1
  • 3
  • 17
2

Maybe I missed the point, but don't you already have your answer? From your post:

order(Mod(x)[1:22050], decreasing=T)[1:10]
[1] 441 881 882 880 883 442 440 879 884 878

Simply collect all values above 500:

junk <- order(Mod(x)[1:22050], decreasing=T)[1:10]
(junk1 <- junk[junk > 500])
[1] 881 882 880 883 879 884 878

To generate the new signal simply repeat what you did to build the original signal:

junk2 <- Sine(0, 1)    
for (i in 1:length(junk1)) {     
    junk2 <- junk2 + Sine(junk1[i], 1)    
}    
junk2.s <- as.vector(junk2$sound)

To keep the values below 500:

(junk3 <- junk[junk <= 500])
[1] 441 442 440
bill_080
  • 4,692
  • 1
  • 23
  • 30
  • Ops, it was too easy to be real XD Thanks for the answer! Just a quick question, that will probably be my next official question here: the resulting sound is *terrible*, not exactly what I was looking for. Do you know a way to improve the fft analysis? Is there a better approach to extract frequencies? – Tommaso Sep 21 '11 at 19:20
  • @Tommaso; I think the "bad" sound is due to the multiple frequencies. From your program, try the following: `play(s1)` `play(s2)` and `play(s3)`. It is the mix of frequencies that causes the "bad" sound. Instead of extracting a range of frequencies, maybe you could pick the middle/median frequency in a range. – bill_080 Sep 21 '11 at 19:50
  • @Tommaso; Oops, ran out of time.... Choosing a median frequency can be done by `(junk1 <- median(junk[junk > 500]))`. – bill_080 Sep 21 '11 at 19:57
  • @Tommaso; After futzing around with `play()`, I noticed some screwy results. To make a long story short, try `play(s3)` and `play(s3/2)`. `s3` is made up of two frequencies. If you build up Z frequencies, divide by Z to play it. – bill_080 Sep 21 '11 at 20:11
  • I came to your same conclusion. The problem is to scale the output in the range `(-1, 1)`. This can be done with the `normalize` function of the library `sound`. I tried to rebuild the original sound, after the fft analysis. Plotting `plot(normalize(s3[1:600]))`, and `plot(normalize(junk4[1:600]+junk2[1:600]))` (junk2 contains frequencies > 500, and junk4 frequencies < 500), will show you a pretty good result (although I have to find a better approximation). Thanks for your help bell! – Tommaso Sep 21 '11 at 20:34
1

look at the 'signal' package on cran, one of the filter functions there should do

Stan
  • 11
  • 1