FFT stm32f407 лишние гармоники в выходном массиве Опции

[yanvasiij]

Доброго времени суток!

Процессор stm32f407, библиотека DSP от Stm32CubeMx. Делаю выборку из 1024 отсчетов встроенного АЦП (12 бит) с частотой 1024 Гц. Когда набирается 1024 отсчета отправляю данные в функции быстрого преобразования Фурье из либы DSP.

Код в более приятном виде можно посмотреть ЗДЕСЬ.

CODE

arm_cfft_radix2_instance_q15 fftInstance;
static u32 curInputCursor=0;
static u32 dataIsReady = 0; /**< Флаг для начала преобразования фурье */
static q15_t input [2048]; /**< входной массив данных */
static q15_t output [1024]; /**< выходной массив данных */

/** @brief прерывание по таймеру с частотой 1024 Гц */
extern "C" void timerInterrupt (void)
{
input[curInputCursor++]=uhADCxConvertedValue; /**< Набираю 1024 замера */
if (curInputCursor>=1024)
{
curInputCursor = 0;
dataIsReady = 1; /**< Флаг, сигнализирующий, что выборка готова */
HAL_TIM_Base_DeInit(&htim1); /**< Как только набралось отключаю таймер, чтобы остановить забор замеров */
}
invertLed(); /**< Инвертирую ногу, чтобы на осциллографе проверить действительно ли частота дискретизации 1024 Гц (да, именно так) */
}

/** @brief задача, которая занимается преобразованием фурье */
static void fftTask (void *p)
{
q15_t maxValue; /* Max FFT value is stored here */
uint32_t maxIndex; /* Index in Output array where max value is */

while(1)
{
if (dataIsReady)
{
dataIsReady = 0;

arm_cfft_radix2_init_q15 (&fftInstance, 1024, 0, 1);
arm_cfft_radix2_q15(&fftInstance, input);
arm_cmplx_mag_q15 (input, output, 1024);
u32 index;
saveFft();
timInit(); /**< Снова запускаю таймер для подготовки следующих 1024 замеров */
}
vTaskDelay(5000);
}
}

Если я все правильно понимаю то в массиве output[] должна получиться частотная характеристика с шагом 1 Гц. При этом output[0] - это постоянная составляющая, все остальное гармоники.
Подаю синусоидальный сигнал 50 Гц с амплитудой 1 вольт, смещенный на 1 вольт. Ожидаю, что в output[0] будет значение постоянной составляющей, а в output[50] амплитуда 50-ой гармоники.
Однако в результате output[0] всегда равен нулю, в output[1] некоторое значение очень похожее на постоянную составлющую, в output[3] тоже непонятно какое значение (181), а в output[50] как и ожидалось значение гармоники (оно смещается если изменить частоту синусоидального сигнала). Вдобавок нет нет да и появятся непонятные значения в разных гармониках (то output[200], то output[511] и т.п. совершенно случайным образом). Добавил RC-цепочку на вход - не помогло.

Буду очень признателен, если кто растолкует, что я делаю неправильно.

Вот здесь картинка с характеристикой (строю на графике массив output[])

[yanvasiij]

Спасибо, теперь всё гораздо понятнее.
А входной сигнал выложите, пожалуйста.

И Вы пишете, что подали 500мВ, а какое напряжение полной шкалы? АЦП двухполярный? Т.е. 1 бит под знак, а значения в доп.коде на выходе? И вообще, лучше давайте продолжать работать не с живыми данными АЦП. Сначала надо отладить всё на идеальной модели входного сигнала, где всё точно.

С какой разрядностью на входе и выходе работает функция arm_cmplx_mag_q15()?

И ещё Вы пишете:
Если не будет выхода за предела разрядной сетки, сдвиг можно сделать хоть до, хоть после функции. Но Вы ещё не назвали разрядность, поэтому как будет известно - будем думать.

И сдвиг вправо, а не влево. Это ошибка в Ваших расчётах или просто опечатка? Если ошибка, то прошу все вышеприведённые цифры повторить заново после исправления.

Функция arm_cmplx_mag_q15() принимает на вход массив, состоящий 16 битных знаковых значений. В описании к функции написано следующее:

The function implements 1.15 by 1.15 multiplications and finally output is converted into 2.14 format.

Я еще так и до конца и не разобрался, что значит форматы 1.15 и 2.14. Но из описания следует, что на вход подаем в формате 1.15, а на выходе имеем 2.14. Более того, входные значения у этой функции в комплексной форме, как я уже и писал - input[n] - действительная часть, input[n+1] - комплексная (n=0,2,4,6,8...). В эту функцию я отправляю массив только после БПФ, предварительно сдвинув каждый член на 10 бит ВЛЕВО, ошибки тут нет:

Код

for (u32 i=0; i<2048; i++) input[i] = input[i]<<10; //Насколько я понимаю -  это сдвиг ВЛЕВО, или я чего то конкретно упустил???

Сделал все сначала на сгенерированном массиве. По-порядку, что я сделал:

1) Сгенерировал 1024 точки синуса с частотой 50 Гц и амплитудой 600. Сигнал двуполярный один бит знак, остальное в дополнительном коде. Вот по такой схеме (с матлабом связываться не стал, проще было прям внутри моей программы его генерировать):

Код

    u32 inputShift=0;
    for (u32 i = 0; i < FFT_LEN; i++)
    {
        float tmp;
        tmp = TEST_AMPL * sin( (TEST_FREQUENCY*2*PI/(FFT_LEN-1)) * i);
        int tmpInt;
        tmpInt = tmp;
        input[inputShift++] = tmpInt;
        input[inputShift++] = 0;
    }

Получилась вот такая картинка:



2) Далее полученные массив отправил на обработку в DSP либу, получилась вот такая характеритика:

[Krys]

Функция arm_cmplx_mag_q15() массив, состоящий 16 битных знаковых значений. В описании к функции написано следующее:

Мы говорим про функцию, которая рассчитывает амплитуду частотных составляющих по их мнимой и действительной частям?

Я еще так и до конца и не разобрался, что значит форматы 1.15 и 2.14. Но из описания следует, что на вход подаем в формате 1.15, а на выходе имеем 2.14.

Ну я Вам про форматы пытался выше объяснить. Если что - спрашивайте конкретные вопросы, поясню более детально.
В рамках конкретной функции такие форматы можно объяснить: если каждую квадратуру выхода БПФ считать пронормированным к 1, то амплитуда в корень из 2 будет больше 1, поэтому при сохранении входного формата возникнет переполнение, следовательно положение точки нужно сдвинуть.
Но на самом деле этого не требуется конкретно для БПФ, т.к. заранее известно, что комплексные числа - это вращающиеся вектора, которые не могут выходить за единичную окружность, тогда и модуль не может превышать 1.

предварительно сдвинув каждый член на 10 бит ВЛЕВО, ошибки тут нет:

А надо-то вправо. БПФ усилил на гейн. Надо вернуть обратно.

1) Сгенерировал 1024 точки синуса с частотой 50 Гц и амплитудой 600

А частота дискретизации тогда какая?

Получилась вот такая картинка:
2) Далее полученные массив отправил на обработку в DSP либу, получилась вот такая характеритика:

По этим картинкам видно, что частота сигнала у Вас не кратна частоте дискретизации. Подайте кратную частоту, так гораздо удобнее сначала разбираться, почему БПФ не дал ожидаемую амплитуду.

[Krys]

А надо-то вправо. БПФ усилил на гейн. Надо вернуть обратно.

Я извиняюсь, запутался. Действительно, здесь надо влево. БПФ усилил, но, чтобы влезть в разрядность сдвинул вправо с потерей младших битов. А чтобы получить число в его реальном масштабе, нужно обратно сдвинуть влево.