Набросаем на плюсах как выглядит комплексное число. Найдем его модуль и аргумент, используя стандартные библиотеки:
Exported from Notepad++
Вывод: 5, 53.1301
Выглядит сложно. Поэтому придумали другой способ вычисления аргумента и, заодно, модуля комплексного числа, названный CORDIC. Это итерационный метод, хорошо описанный в wiki. Я попробую объяснить его немного с другой стороны.
Пусть у нас есть комплексное число Z = 3+4j. На комплексной плоскости оно выглядит так:
Аргумент этого числа рассчитали в программе выше и он равен 53.1301, или для простоты 53 градуса. Для этого вектора выполняются следующие равенства:
Re = |Z|*cos(53) = 3
Im = |Z|*sin(53) = 4, где |Z| - модуль комплексного числа.
Повернем этот вектор на L градусов, тогда:
Re_new = |Z|*cos(53+L) = |Z|*(cos(53)*cos(L) - sin(53)*sin(L))
Im_new = |Z|*sin(53+L) = |Z|*(sin(53)*cos(L) + cos(53)*sin(L)).
Тут надо вспомнить школу и тригонометрические формулы. Just google it. Если раскрыть скобки и воспользоваться подстановкой |Z|*cos(53) = Re и |Z|*sin(53) = Im, то получится:
Re_new = Re*cos(L) - Im*sin(L)
Im_new = Im*cos(L) + Re*sin(L)
Получается логичная вещь: значения поворачиваемого вектора зависят от первоначального положения вектора и угла, на который вращают вектор. При этом нет разницы повернем мы вектор один раз на L градусов, или два раза на L/2 градусов или N раз на L/N градусов.Если мы вращаем на отрицательный угол, то поменяется лишь знак:
Re_new = Re*cos(-L) - Im*sin(-L) = Re*cos(L) + Im*sin(L)
Im_new = Im*cos(-L) + Re*sin(-L) = Im*cos(L) - Re*sin(L)
Теперь рассмотрим такую ситуацию: у нас на входе комплексное число с произвольными Re и Im. Как определить аргумент такого числа? Мы просто будем поворачивать этот вектор до тех пор, пока Re или Im не станет равной нулю. При повороте вектора его модуль не меняется и в случае, когда одна из составляющих комплексного числа равно нулю, то другая равна его модулю. Вот и весь смысл метода CORDIC. Теперь рассмотрим более детально реализацию.
Видно, что для расчета Re_new и Im_new в ПЛИС нужно хранить две таблицы: синусов и косинусов. Их количество можно уменьшить вдвое, если поделить все на cos(L):
Re_new/cos(L) = Re - Im*tan(L)
Im_new/cos(L) = Im + Re*tan(L)
Видно, что теперь нужно хранить таблицу только из тангенсов углов. Перепишем выражение выше таким образом:
Re_new_сordic = Re - Im*tan(L)
Im_new_cordic = Im + Re*tan(L)
Проблему с делением на cos(L) на каждом ходу решается довольно просто: если количество итераций(поворотов вектора) сделать константой, то в результате N итераций соотношение между Re_new и Re_new_cordic будет таким:
Re_new = Re_new_cordic * cordic_const, где cordic_const равна произведению косинусов углов поворота и является константой. То есть мы спокойно считаем N итераций результат Re_new_cordic, а потом конечный результат просто умножаем на константу. Стоит отметить, что cordic_const принято называть cordic gain.
Ну и последнее, что важно отметить: на какие углы вращать вектор? Представим, что аргумент равновероятно лежит в диапазоне 0..90 градусов. Логично воспользоваться бинарным поиском и делить отрезки все время пополам, тогда углы будут 45, 22, 11, 5, 2, 1. В этом случае cordic_const = cos(45)*cos(22)*cos(11)*cos(5)*cos(2)*cos(1) и таблица тангенсов состоит из шести значений. Но есть и другой вариант, который подходит больше. Нативным для Verilog являются целочисленные переменные. А для целочисленных переменных деление на степень двойки равносильно операции сдвига. Поэтому нужно подобрать такие углы, тангенс которых был бы равен 1/(2^k), k = 0..N-1, где N - количество итераций. Если все рассуждения перенести в код, то получится такая программа:
Exported from Notepad++
Вывод в файл console.txt:
tangens angles:
45
26.5651
14.0362
7.12502
3.57633
1.78991
0.895174
0.447614
CORDIC gain:
0.607259
Re of Z: 3
Im of Z: 4
CORDIC arg: 57.1128 real: 53.1301
CORDIC module: 5 real: 5
#include <iostream>
#include <math.h>
#include <limits.h>
using namespace std;
const double PI = M_PI;
class complex {
private:
int re;
int im;
public:
complex(){re = 0; im = 0;}
complex(int r, int i){re = r; im = i;}
int get_re(){return re;}
int get_im(){return im;}
void set_re(int r) {re = r;}
void set_im(int i) {im = i;}
int get_abs_val(){return sqrt(re*re+im*im);}
float get_arg();
};
float complex::get_arg() {
float f_re = float(re);
float f_im = float(im);
float arg = atan2f(f_im,f_re)*180/M_PI;
return arg;
}
int main() {
complex *a = new complex(3,4);
cout << a->get_abs_val() << endl;
cout << a->get_arg() << endl;
return 0;
}
Вывод: 5, 53.1301
Как перевести на Verilog вычисление модуля уже разобрались тут и тут. Остается вопрос с фазой. Определение из wiki выглядит следующим образом:
Пусть у нас есть комплексное число Z = 3+4j. На комплексной плоскости оно выглядит так:
Re = |Z|*cos(53) = 3
Im = |Z|*sin(53) = 4, где |Z| - модуль комплексного числа.
Повернем этот вектор на L градусов, тогда:
Re_new = |Z|*cos(53+L) = |Z|*(cos(53)*cos(L) - sin(53)*sin(L))
Im_new = |Z|*sin(53+L) = |Z|*(sin(53)*cos(L) + cos(53)*sin(L)).
Тут надо вспомнить школу и тригонометрические формулы. Just google it. Если раскрыть скобки и воспользоваться подстановкой |Z|*cos(53) = Re и |Z|*sin(53) = Im, то получится:
Re_new = Re*cos(L) - Im*sin(L)
Im_new = Im*cos(L) + Re*sin(L)
Получается логичная вещь: значения поворачиваемого вектора зависят от первоначального положения вектора и угла, на который вращают вектор. При этом нет разницы повернем мы вектор один раз на L градусов, или два раза на L/2 градусов или N раз на L/N градусов.Если мы вращаем на отрицательный угол, то поменяется лишь знак:
Re_new = Re*cos(-L) - Im*sin(-L) = Re*cos(L) + Im*sin(L)
Im_new = Im*cos(-L) + Re*sin(-L) = Im*cos(L) - Re*sin(L)
Теперь рассмотрим такую ситуацию: у нас на входе комплексное число с произвольными Re и Im. Как определить аргумент такого числа? Мы просто будем поворачивать этот вектор до тех пор, пока Re или Im не станет равной нулю. При повороте вектора его модуль не меняется и в случае, когда одна из составляющих комплексного числа равно нулю, то другая равна его модулю. Вот и весь смысл метода CORDIC. Теперь рассмотрим более детально реализацию.
Видно, что для расчета Re_new и Im_new в ПЛИС нужно хранить две таблицы: синусов и косинусов. Их количество можно уменьшить вдвое, если поделить все на cos(L):
Re_new/cos(L) = Re - Im*tan(L)
Im_new/cos(L) = Im + Re*tan(L)
Видно, что теперь нужно хранить таблицу только из тангенсов углов. Перепишем выражение выше таким образом:
Re_new_сordic = Re - Im*tan(L)
Im_new_cordic = Im + Re*tan(L)
Проблему с делением на cos(L) на каждом ходу решается довольно просто: если количество итераций(поворотов вектора) сделать константой, то в результате N итераций соотношение между Re_new и Re_new_cordic будет таким:
Re_new = Re_new_cordic * cordic_const, где cordic_const равна произведению косинусов углов поворота и является константой. То есть мы спокойно считаем N итераций результат Re_new_cordic, а потом конечный результат просто умножаем на константу. Стоит отметить, что cordic_const принято называть cordic gain.
Ну и последнее, что важно отметить: на какие углы вращать вектор? Представим, что аргумент равновероятно лежит в диапазоне 0..90 градусов. Логично воспользоваться бинарным поиском и делить отрезки все время пополам, тогда углы будут 45, 22, 11, 5, 2, 1. В этом случае cordic_const = cos(45)*cos(22)*cos(11)*cos(5)*cos(2)*cos(1) и таблица тангенсов состоит из шести значений. Но есть и другой вариант, который подходит больше. Нативным для Verilog являются целочисленные переменные. А для целочисленных переменных деление на степень двойки равносильно операции сдвига. Поэтому нужно подобрать такие углы, тангенс которых был бы равен 1/(2^k), k = 0..N-1, где N - количество итераций. Если все рассуждения перенести в код, то получится такая программа:
#include <iostream>
#include <fstream>
#include <math.h>
#include <limits.h>
using namespace std;
const double PI = M_PI;
#define N 8 //number steps
class complex {
private:
int re;
int im;
public:
complex(){re = 0; im = 0;}
complex(int r, int i){re = r; im = i;}
int get_re(){return re;}
int get_im(){return im;}
void set_re(int r) {re = r;}
void set_im(int i) {im = i;}
int get_abs_val(){return sqrt(re*re+im*im);}
float get_arg();
};
float complex::get_arg() {
float f_re = float(re);
float f_im = float(im);
float arg = atan2f(f_im,f_re)*180/M_PI;
return arg;
}
int main() {
ofstream console("console.txt", ios::out);
float d_angle[N];
float sum_angles = .0f;
console << "tangens angles:\n";
for(int i = 0; i < N; i++) {
d_angle[i] = atan(1.0/pow(2,i)) * (180/PI);
sum_angles += d_angle[i];
console << d_angle[i] << endl;
}
float cordic_gain = .0;
for(int i =0; i < N; i++) {
if(i == 0) {
cordic_gain = cos(d_angle[i]*PI/180.0);
}
else {
cordic_gain *= cos(d_angle[i]*PI/180.0);
}
}
console << "CORDIC gain:\n" << cordic_gain << endl;
short Re[N+1];
short Im[N+1];
float arg = .0f;
complex *a = new complex(3,4);
Re[0] = a->get_re();
console << "Re of Z: " << Re[0] << endl;
Im[0] = a->get_im();
console << "Im of Z: " << Im[0] << endl;
for(int i = 0; i < N; i++) {
if(Im[i] >= 0) {
Re[i+1] = Re[i] + (Im[i] >> i);
Im[i+1] = Im[i] - (Re[i] >> i);
arg -= d_angle[i];
}
else {
Re[i+1] = Re[i] - (Im[i] >> i);
Im[i+1] = Im[i] + (Re[i] >> i);
arg += d_angle[i];
}
cout << i << "\t" << Re[i+1] << "\t" << Im[i+1] << "\t" << arg << endl;
}
cout << "\nCORDIC:\t" << -arg << "\t real: " << a->get_arg() << endl;
console << "\nCORDIC arg:\t" << -arg << "\t real: " << a->get_arg() << endl;
int my_abs = Re[N] * cordic_gain;
cout << "\nCORDIC:\t" << my_abs << "\t real: " << a->get_abs_val() << endl;
console << "\nCORDIC module:\t" << my_abs << "\t real: " << a->get_abs_val() << endl;
console.close();
return 0;
}
Вывод в файл console.txt:
tangens angles:
45
26.5651
14.0362
7.12502
3.57633
1.78991
0.895174
0.447614
CORDIC gain:
0.607259
Re of Z: 3
Im of Z: 4
CORDIC arg: 57.1128 real: 53.1301
CORDIC module: 5 real: 5
Далее будем разбираться с результатами и переносить код на Verilog
Комментариев нет:
Отправить комментарий