Как получить частоту гармоники спектра?

123
08 июля 2019, 06:20

Генерирую синусоиду

for(quint32 i = 0; i < ALen; ++i) ABuffer[i] = AAmplitude * qSin(i * M_2_PI * ASigFreq/ADiscrFreg);

С помощью библиотеки FFTW получаю его спектр в комплексном виде

// подготовка плана преобразования
fftw_plan myPlan = fftw_plan_dft_r2c_1d(static_cast<int>(SampleCount), z, y, FFTW_ESTIMATE);
// расчет ДПФ
fftw_execute(myPlan);

где SampleCount - количество отсчетов

z - массив описывающий синусоиду

y - массив результатов

Получаю массив амплитуд в вещественном виде

for (int i = 0; i < fftSampleCount; ++i) fftSamples[i] = sqrt(y[i][0] * y[i][0] + y[i][1] * y[i][1]);

Вычисляю "пик" и его положение в спектре

auto GarmAmplMax = std::max_element(fftSamples,(fftSamples+ADataLen));
auto GarmNum = std::distance(fftSamples,GarmAmplMax);

GarmAmplMax - амплитуда пика

GarmNum - позиция в спектре fftSamples

Внимание! Вопрос: как рассчитать частоту для полученного пика?

Код функции целиком

void MainWindow::on_actGetData_triggered()
{
ui->pbRecord->setEnabled(false);
// производим запись, а не генерацию звука
qint64 SampleCount;
double * x = nullptr;
double * z = nullptr;
QString WindowName;
auto getCoeff = [&](int i ) -> double
{
    switch (FWindowIndex)
    {
        case 0 : // без окна
            WindowName = "Без окна";
            return 1;
        case 1 : // Hann
            WindowName = "Окно Ханна";
            return 0.5*(1 - qCos((2 * M_PI * i)/(SampleCount-1)));
        case 2 : // Gausse
            {
                WindowName = "Окно Гауса";
                double a = (SampleCount-1)/2;
                double t = (i - a)/(0.5 * a);
                return qExp(-(t*t)/2);
            };
        case 3 : // Hamming
            WindowName = "Окно Хэммнга";
            return 0.54 - 0.46 * qCos((M_2_PI * i)/(SampleCount - 1));
    }
    return 0;
};
int devIndex = ui->cbDevice->currentIndex();
if(ui->cbDevice->currentIndex())
{
    QAudioFormat format;
    format.setSampleRate(ui->spbDescrFreq->value());
    format.setChannelCount(1);
    format.setSampleSize(16);
    format.setCodec("audio/pcm");
    format.setByteOrder(QAudioFormat::LittleEndian);
    format.setSampleType(QAudioFormat::SignedInt);
    QList<QAudioDeviceInfo> list = QAudioDeviceInfo::availableDevices(QAudio::AudioInput);
    QAudioDeviceInfo info = list.at(devIndex - 1);
    if(!info.isFormatSupported(format))
    {
        addMessage("Заданный формат записи не поддерживается");
        ui->pbRecord->setEnabled(true);
        return;
    };
    FRecorder = new QAudioInput(info,format, this);
    FBuffer->open(QIODevice::ReadWrite);
    FRecorder->start(FBuffer);
    addMessage("Запись запущена.");
    QTest::qWait(ui->spbDuration->value());
    addMessage("Запись окончена. Начата обработка данных.");
    SampleCount = (*FBuffer).size()/2;
    qint16 * SamplesList = reinterpret_cast<qint16 *>((*FBuffer).buffer().data());
    addMessage(QString("Сделано %1 отсчетов").arg(SampleCount));
    x = new double [SampleCount];
    z = new double [SampleCount];
    // накладываем окно
    for (int i = 0; i < SampleCount; ++i)
    {
        x[i] = SamplesList[i];
        z[i] = getCoeff(i) * x[i];
    }
    FBuffer->buffer().clear();
    delete FRecorder;
    FRecorder = nullptr;
}
else // генерация звука
{                       // секунды                      частота
    SampleCount = (ui->spbDuration->value()/1000) * ui->spbDescrFreq->value();
    addMessage("Начата генерация данных.");
    x = new double [SampleCount];
    z = new double [SampleCount];
    generateSinWav(x,                          // приемник данных
                   SampleCount,                // количество отсчетов
                   ui->spbDescrFreq->value(),  // частота дискретизации
                   ui->spdSignalFreq->value(), // частота синусоиды
                   20000);                     // амплитуда сигнала
    for (int i = 0; i < SampleCount; ++i)
    {
        z[i] = getCoeff(i) * x[i];
        //addMessage(QString("x[%1] = %2").arg(i).arg(x[i]));
    }
    addMessage("Данные сгенерированы. Начата их обработка");
}
qint64 fftSampleCount = SampleCount / 2 + 1;
fftw_complex * y = new fftw_complex[fftSampleCount];
// подготовка плана преобразования
fftw_plan myPlan = fftw_plan_dft_r2c_1d(static_cast<int>(SampleCount), z, y, FFTW_ESTIMATE);
// расчет ДПФ
fftw_execute(myPlan);
fftw_destroy_plan(myPlan);
double * fftSamples = new double[fftSampleCount];
for (int i = 0; i < fftSampleCount; ++i) fftSamples[i] = sqrt(y[i][0] * y[i][0] + y[i][1] * y[i][1]); //получаем модуль комплексного числа - амплитуду сигнала
addMessage("обработка закончена");
ChartForm * Form = new ChartForm;
Form->setAttribute(Qt::WA_DeleteOnClose, true);
Form->showSourceSignalChart(x,SampleCount,ui->spdSignalFreq->value());
Form->showSignalChart(z,SampleCount,WindowName);
Form->showSpectrChart(fftSamples,fftSampleCount,ui->spbDescrFreq->value());
Form->show();
if(fftSamples) delete [] fftSamples;
if(z) delete [] z;
if(y) delete [] y;
if(x) delete [] x;
ui->pbRecord->setEnabled(true);
}
Answer 1

Отсчёт с номером N/2 (где N - SampleCount), соответствует половинной частоте дискретизации исходного сигнала в герцах.

Например, для звука 44.1 кГц и интервала преобразования размером 4096:
22050 Гц соответствует Result[2048]
нулевой отсчёт соответствует нулевой частоте (постоянная составляющая)
первый отсчёт - 10.8 Гц

Таким образом, частота, соответствующая i-му отсчёту результата:

 f = Freq_Discr * i / N
READ ALSO
Присваивание массива

Присваивание массива

Вопрос в следующем: если присваивать в цикле:

132
создание сервиса с помощью Java SPI

создание сервиса с помощью Java SPI

Необходимо сделать расширяемое за счет плагинов приложениеДля создания подобной системы проделал следующие шаги:

109
Как изменит background правильно?

Как изменит background правильно?

За ранее я извиняюсь, Я не силён на русском языкеЯ создаю приложения

131
Менять фон layout&#39;а в виджет программно

Менять фон layout'а в виджет программно

Всем приветВозникла такая проблема, нужно при нажатии кнопки в самом приложении менять фон layout'а в виджете

119