Pitch Estimation Algorithm

Started by
0 comments, last by alvaro 12 years ago
Hi,

I'm trying to implement an algorithm that, given human singing samples captured from a mic, tell the singing note.
I'm no math expert, but I've read quite a lot of material on this topic the past weeks. I have an algorithm working but with a few problems, that I don't know how to solve them.

Basically what I do is:

1) Apply a hamming window to raw input data. (1024 samples from mono mic input at 44100Hz rate)
2) Apply FFT
3) Get magnitude and true frequency (frequency computed from the bin frequency + phase offset) for each bin
4) Return the frequency of the bin with greatest magnitude.

I'm using this tool to test my results: http://www.seventhst...tuningfork.html
I place the mic next to the speaker to capture input. I've noticed that this estimation depends on how far is the mic from my speakers. If I put it right in front of it I get right results for tones greater than D2 (~293Hz) even with speakers volumes at the minimum. Below that freq it gives me completely wrong values. If I move the mic 5 inches away of the speakers I start getting wrong values from below G2 (~392Hz) even if the speaker volumes are at the maximum.

It seems something is either wrong with my algorithm or my mic, or both.

The algorithm follows, pehaps you could shed some light on it:


// Using portaudio to capture mic Input. This callback is called each frame with inputBuffer filled with raw data.
int paCallback(const void* inputBuffer, void* outputBuffer, unsigned long framesPerBuffer,
const PaStreamCallbackTimeInfo* timeInfo, PaStreamCallbackFlags statusFlags,
void* userData)
{
float* input = (float*) inputBuffer;
double* data = (double*)userData;
if(input != NULL)
{
for(unsigned int i = 0; i < FFT_N; i++)
{
// apply hamming window
data[2 * i] = input * fc.getWindow(i);
data[2 * i + 1] = 0.0;
}
double freq = 0.0;
double ampl = 0.0;

// fc is a global instance of Analyzer class
fc.analyze(data, FFT_N, SAMPLE_RATE, 1, &freq, &ampl);
double db = log10(ampl);
printf("%-10s | %8.3lfHz | %5.2lfdB | %lf\n", Note::noteName(freq), freq, db, ampl);
}

return 0;
}


void FrequencyCounter::analyze(double* data, unsigned long nn, double sampleRate, int overlapFactor, double* outFreq, double* outMagnitude)
{
// daniel-lanzos algorithm (reverse binary reindexing dark magic)
fft(data, FFT_N);

// Precalculated constants
const double freqPerBin = sampleRate / FFT_N;
const double stepSize = FFT_N / overlapFactor;
const double expectPhaseDiff = 2.0 * M_PI * stepSize / FFT_N;
double real = 0.0;
double imag = 0.0;
double phase = 0.0;
double delta = 0.0;
long qpd = 0;

const size_t iMax = std::min(size_t(FFT_N / 2), size_t(FFT_MAXFREQ / freqPerBin));
for (size_t i = 0; i < iMax; ++i)
{
real = data[2*i];
imag = data[2*i+1];

phase = atan2(imag, real);

// process phase difference
delta = phase - mFFTLastPhase;
mFFTLastPhase = phase;

// subtract expected phase difference
delta -= i * expectPhaseDiff;

/* map delta phase into +/- Pi interval */
qpd = delta / M_PI;
if (qpd >= 0)
{
qpd += qpd & 1;
}
else
{
qpd -= qpd & 1;
}
delta -= M_PI * static_cast<double>(qpd);

/* get deviation from bin frequency from the +/- Pi interval */
delta = overlapFactor * delta / (2.0 * M_PI);

// true frequency
data[2 * i] = (i + delta) * freqPerBin;

// magnitude
data[2 * i + 1] = 2.0 * sqrt(real * real + imag * imag);
}

unsigned int maxI = 0;
double maxMag = data[1];
for(unsigned int i = 0; i < iMax; i++)
{
if(data[2 * i + 1] > maxMag)
{
maxI = i;
maxMag = data[2 * i + 1];
}
}

if(outFreq != NULL)
{
*outFreq = data[2 * maxI];
}
if(outMagnitude != NULL)
{
*outMagnitude = data[2 * maxI + 1];
}
}


Thanks in advance.
Advertisement
The bin with the largest magnitude is not necessarily the fundamental frequency. You probably need a better algorithm. Perhaps you can start here:
http://en.wikipedia.org/wiki/Pitch_detection_algorithm
http://en.wikipedia.org/wiki/Frequency_estimation

It's probably a good idea to make some pictures of the spectrum with real input, and then try several algorithms on that data.

This topic is closed to new replies.

Advertisement