Не получение точных данных при обратном БПФ

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

Код:

const int fft_siz = 512;
const int inverse = 1;

kiss_fft_cpx* in = (kiss_fft_cpx*)malloc(sizeof(kiss_fft_cpx) * fft_siz);
kiss_fft_cpx* out = (kiss_fft_cpx*)malloc(sizeof(kiss_fft_cpx) * fft_siz);
kiss_fft_cpx* rec = (kiss_fft_cpx*)malloc(sizeof(kiss_fft_cpx) * fft_siz);
kiss_fft_cfg cfg = kiss_fft_alloc(fft_siz, !inverse, NULL, NULL);
kiss_fft_cfg icfg = kiss_fft_alloc(fft_siz, inverse, NULL, NULL);

srand((unsigned int)time(NULL));
for(int i = 0; i < fft_siz; i++)
{
    in[i].r = rand() % 256;
    in[i].i = rand() % 256;
}

kiss_fft(cfg, in, out);

// scaling
for(int i = 0; i < fft_siz; i++)
{
    out[i].r /= fft_siz;
    out[i].i /= fft_siz;
}

kiss_fft(icfg, out, rec);

unsigned int count = 0;
for(int i = 0; i < fft_siz; i++)
    if(in[i].r != rec[i].r)
    {
        count++;
        printf( "in[%3d].r does not match rec[%3d].r :: %3d :: %f\n",
                i, i, count, in[i].r - rec[i].r);
    }
    else if(in[i].i != rec[i].i)
    {
        count++;
        printf( "in[%3d].i does not match rec[%3d].i :: %3d :: %f\n",
                i, i, count, in[i].i - rec[i].i);
    }

free(in);
free(out);
free(rec);
free(cfg);
free(icfg);

kiss_fft_cleanup();

Выход:

in[  0]:     71.000000       85.000000 -- out[  0]:    127.095703      124.541016
in[  1]:    248.000000       27.000000 -- out[  1]:     -7.083314        0.072701
in[  2]:     64.000000       18.000000 -- out[  2]:     -3.770610        2.682554
in[  3]:      6.000000       96.000000 -- out[  3]:     -7.929140       -2.897723
in[  4]:     98.000000       23.000000 -- out[  4]:     -0.719621       -5.854260
in[  5]:    250.000000      188.000000 -- out[  5]:      0.397226       -1.248124
in[  6]:    231.000000        3.000000 -- out[  6]:     -7.934285       -2.367196
in[  7]:      6.000000      105.000000 -- out[  7]:     -0.317480       -2.955601
in[  8]:    172.000000      143.000000 -- out[  8]:     -4.236186        3.911616
in[  9]:     16.000000      134.000000 -- out[  9]:     -0.162577       -5.353521
in[ 10]:    230.000000      112.000000 -- out[ 10]:     -4.703711        7.791993
in[ 11]:      5.000000       26.000000 -- out[ 11]:     -2.636305        0.188381
in[ 12]:     16.000000      127.000000 -- out[ 12]:      1.137413        4.576081
in[ 13]:    112.000000       86.000000 -- out[ 13]:      0.978051       -0.408992
in[ 14]:     40.000000       23.000000 -- out[ 14]:      5.231920       -2.347566
in[ 15]:     75.000000       26.000000 -- out[ 15]:      0.009981       -2.091559

note                                ::count::difference
--------------------------------------------------------
in[  1].r does not match rec[  1].r ::   1 :: -0.000031
in[  3].r does not match rec[  3].r ::   2 :: -0.000015
in[  4].i does not match rec[  4].i ::   3 :: -0.000004
in[  6].i does not match rec[  6].i ::   4 :: -0.000008
in[  7].r does not match rec[  7].r ::   5 :: -0.000002
in[  9].r does not match rec[  9].r ::   6 :: -0.000015
in[ 11].r does not match rec[ 11].r ::   7 :: -0.000015
in[ 12].r does not match rec[ 12].r ::   8 :: -0.000015
in[ 13].i does not match rec[ 13].i ::   9 :: -0.000008
in[ 14].i does not match rec[ 14].i ::  10 :: 0.000008
in[ 15].r does not match rec[ 15].r ::  11 :: -0.000015

Отладка. Если вы пойдете вниз, вы увидите, что там 317 несоответствий. Я также вывожу разницу между значениями, т.е. (in[].r - rec[].r) или (in[].i - rec[].i).

Далее я показываю входные данные, где белые точки представляют реальную часть, а красные точки — мнимую часть.

введите здесь описание изображения

Это выходные данные БПФ, представленные фиолетовым, а реконструированные данные — белым и красным.

введите здесь описание изображения

введите здесь описание изображения

Заметили небольшую разницу? Я предполагаю, что это связано с точностью с плавающей запятой. Как я могу решить эту проблему, чтобы получить точно такие же входные данные, на которых я использовал БПФ?

Редактировать:

I noticed that the range off error in my case is something like ]0 , 0.002]. So as a workaround, I rounded the reconstructed data and got a good result. But still... This only works if the fractional part of my numbers is 0.0.


person Jonas    schedule 30.06.2013    source источник


Ответы (1)


Вы правильно догадываетесь; числа с плавающей запятой имеют ограниченную точность, и часть входной точности будет потеряна в процессе. Если вам нужен точный ввод, вам понадобится какая-нибудь библиотека произвольной точности с плавающей запятой (например, библиотека GNU MP< /а>).

person Drew McGowen    schedule 30.06.2013
comment
мда, было бы больно модифицировать библиотеку, используя другую, зная, что я не писал ни одну из них. - person Jonas; 30.06.2013
comment
Библиотеки произвольной точности недостаточно для получения точных результатов. Многие коэффициенты, используемые в БПФ, иррациональны (более того, трансцендентны). Таким образом, они не могут быть точно представлены с какой-либо конечной точностью или базой. Большая точность может только уменьшить ошибку, но не устранить ее. Тем не менее, как правило, нет причин для получения точных результатов в БПФ. Любые данные, полученные из реального мира (например, через микрофон), имеют ошибки измерения, поэтому любая их обработка должна допускать ошибки. Вопрос только в том, насколько, удерживая шум ниже приемлемого уровня. - person Eric Postpischil; 30.06.2013
comment
Значит, потеря некоторой точности неизбежна?! - person Jonas; 30.06.2013
comment
Что ж, если вам абсолютно нужны точные результаты, вы можете использовать какую-нибудь символическую математическую библиотеку. - person Drew McGowen; 30.06.2013