Инструменты пользователя

Инструменты сайта


algorithm

Алгоритмы

Тригонометрия

CORDIC

Источник

CORDIC - это простой и эффективный алгоритм вычисления синуса и косинуса значения с использованием только базовой арифметики (сложение, вычитание и сдвиги). Ниже представлен очень простой ANSI C код реализующий CORDIC для арифметики с фиксированной точкой. Он основан на определениях из прекрасной книги FXTBook. Прочитайте её, если вам интересны подробности.

Приведённая ниже программа генерирует заголовочный файл с константами и кодом для CORDIC вычислений, используя различные размеры слова. Программа выводит заголовочный файл на стандартный вывод. Следует отметить, что данный код использует математическую библиотеку для расчёта констант, то есть его следует запускать на устройстве оборудованном FPU, а затем скомпилировать получившийся заголовочный файл для вашего целевого устройства. Укажите количество бит, которые должны быть использованы, это приведёт к генерации таблиц, которые соответствуют 1.0 == 2^(bitsize-2). Например, для 16 бит, получившийся заголовочный файл будет использовать представление фиксированной точки вида 2.14. Это определяет максимум точности. Вы можете изменить mul ниже на другое значение, если у вас уже есть необходимое представление фиксированной точки и вы желаете его использовать.

gentable.c
#include <math.h>
#define M_PI 3.1415926535897932384626
#define K1 0.6072529350088812561694
 
int main(int argc, char **argv)
{
    int i;
    int bits = 32; // number of bits 
    int mul = (1<<(bits-2));
 
 
    int n = bits; // number of elements. 
    int c;
 
    printf("//Cordic in %d bit signed fixed point math\n", bits);
    printf("//Function is valid for arguments in range -pi/2 -- pi/2\n");
    printf("//for values pi/2--pi: value = half_pi-(theta-half_pi) and similarly for values -pi---pi/2\n");
    printf("//\n");
    printf("// 1.0 = %d\n", mul);
    printf("// 1/k = 0.6072529350088812561694\n");
    printf("// pi = 3.1415926535897932384626\n");
 
    printf("//Constants\n");
    printf("#define cordic_1K 0x%08X\n", (int)(mul*K1));
    printf("#define half_pi 0x%08X\n", (int)(mul*(M_PI/2)));
    printf("#define MUL %f\n", (double)mul);
    printf("#define CORDIC_NTAB %d\n", n);
 
    printf("int cordic_ctab [] = {");
    for(i=0;i<n;i++)
    {
        c = (atan(pow(2, -i)) * mul);
        printf("0x%08X, ", c);
    }
    printf("};\n\n");
 
    //Print the cordic function
    printf("void cordic(int theta, int *s, int *c, int n)\n{\n  int k, d, tx, ty, tz;\n");
    printf("  int x=cordic_1K,y=0,z=theta;\n  n = (n>CORDIC_NTAB) ? CORDIC_NTAB : n;\n");
    printf("  for (k=0; k<n; ++k)\n  {\n    d = z>>%d;\n", (bits-1));
    printf("    //get sign. for other architectures, you might want to use the more portable version\n");
    printf("    //d = z>=0 ? 0 : -1;\n    tx = x - (((y>>k) ^ d) - d);\n    ty = y + (((x>>k) ^ d) - d);\n");
    printf("    tz = z - ((cordic_ctab[k] ^ d) - d);\n    x = tx; y = ty; z = tz;\n  }  \n *c = x; *s = y;\n}\n");
}

CORDIC функция будет ожидать на входе значение из диапазона -pi/2 .. pi/2. Другие значения можно легко вычислить. Для pi/2–pi результатом будет (pi/2-(a-pi/2)) и аналогично для -pi—pi/2. Любой другой угол может быть отображён на этот диапазон взяв угол по модулю pi.

Предположим, что вы используете 32 бита и сохранили вывод предыдущей программы в файл cordic-32bit.h, вы получите подобный файл:

cordic-32bit.h
//Cordic in 32 bit signed fixed point math
//Function is valid for arguments in range -pi/2 -- pi/2
//for values pi/2--pi: value = half_pi-(theta-half_pi) and similarly for values -pi---pi/2
//
// 1.0 = 1073741824
// 1/k = 0.6072529350088812561694
// pi = 3.1415926535897932384626
//Constants
#define cordic_1K 0x26DD3B6A
#define half_pi 0x6487ED51
#define MUL 1073741824.000000
#define CORDIC_NTAB 32
int cordic_ctab [] = {0x3243F6A8, 0x1DAC6705, 0x0FADBAFC, 0x07F56EA6, 0x03FEAB76, 0x01FFD55B, 
0x00FFFAAA, 0x007FFF55, 0x003FFFEA, 0x001FFFFD, 0x000FFFFF, 0x0007FFFF, 0x0003FFFF, 
0x0001FFFF, 0x0000FFFF, 0x00007FFF, 0x00003FFF, 0x00001FFF, 0x00000FFF, 0x000007FF, 
0x000003FF, 0x000001FF, 0x000000FF, 0x0000007F, 0x0000003F, 0x0000001F, 0x0000000F, 
0x00000008, 0x00000004, 0x00000002, 0x00000001, 0x00000000, };
 
void cordic(int theta, int *s, int *c, int n)
{
  int k, d, tx, ty, tz;
  int x=cordic_1K,y=0,z=theta;
  n = (n>CORDIC_NTAB) ? CORDIC_NTAB : n;
  for (k=0; k<n; ++k)
  {
    d = z>>31;
    //get sign. for other architectures, you might want to use the more portable version
    //d = z>=0 ? 0 : -1;
    tx = x - (((y>>k) ^ d) - d);
    ty = y + (((x>>k) ^ d) - d);
    tz = z - ((cordic_ctab[k] ^ d) - d);
    x = tx; y = ty; z = tz;
  }  
 *c = x; *s = y;
}

Очевидно, что 16 битный CORDIC алгоритм может быть создан также легко (результатом будет cordic-16bit.h). Следующий файл показывает пример использования CORDIC для простой тестовой функции, которая сравнивает результаты со стандартной математической реализацией.

cordic-test.c
#include "cordic-32bit.h"
#include <math.h> // for testing only!
 
//Print out sin(x) vs fp CORDIC sin(x)
int main(int argc, char **argv)
{
    double p;
    int s,c;
    int i;    
    for(i=0;i<50;i++)
    {
        p = (i/50.0)*M_PI/2;        
        //use 32 iterations
        cordic((p*MUL), &s, &c, 32);
        //these values should be nearly equal
        printf("%f : %f\n", s/MUL, sin(p));
    }       
}

CORDIC алгоритм не самый точный, но он прост в реализации на железе с ограниченными ресурсами (например, без аппаратной реализации умножения или в случае малого ПЗУ).

algorithm.txt · Последние изменения: 2017/02/11 11:49 — Ruslan Popov