Пол года назад я наткнулся в сети вот на это видео.
Первой мыслью было то, что это очень круто и у меня такое никогда не получится повторить. Шло время, читались статьи, изучались методы и я искал примеры реализации подобного, но к моему огорчению, в сети ничего конкретного не находилось. Наткнувшись однажды на вычисления тригонометрических функций с использованием алгоритмов CORDIC, я решил попробовать создать свою собственную вращалку изображения на ПЛИС.
CORDIC
Итак, CORDIC — это аббревиатура от COordinate Rotation DIgital Computer.
Это мощный инструмент для вычисления гиперболических и тригонометрических функций. Большинство алгоритмов CORDIC работают методом последовательного приближения и не очень сложны в реализации как на языках программирования высокого уровня, так и на HDL. Я не стану заострять внимание на математике метода, читатель может ознакомиться с ним в сети или по ссылкам ниже.
В свободном доступе мне попалась вот эта реализация алгоритма CORDIC на языке verilog. Данное ядро работает в 2-х режимах: Rotate и Vector. Для наших целей подходит режим Rotate. Он позволяет вычислять значения функций sin и cos от заданного угла в радианах или градусах. Библиотеку можно сконфигурить как в конвейерном, так и в комбинационном варианте. Для наших целей подходит конвейер, у него самое большое Fmax. Он выдаст значения синуса и косинуса с задержкой в 16 тактов.
В RTL Viewer-e модуль CORDIC отображается состоящим из 16 однотипных блоков:
Каждый из которых принимает на вход данные с предыдущего и выходами подключен ко входам следующего. Выглядит он так:
Ядро библиотеки работает только в первом квадранте, а это значит что оставшиеся три нам придётся вычислять самим вычитая pi/2 и меняя знак.
Выбранный мной подход не является очень правильным т.к. качество вращаемого изображения оставляет желать лучшего. Это происходит по причине расчета координат на-лету, без применения дополнительной буферизации данных и последовательного вычисления координат за несколько проходов, как это делается в Shear.
Первой инстанцией нашего вращателя является блок расчёта квадранта и угла поворота. Угол поворота инкрементируется каждый новый кадр на 1 градус. По достижению угла 90 градусов, квадрант меняется на следующий по очереди, а угол либо сбрасывается в ноль, либо декрементируется на 1 градус каждый новый кадр.
Выглядит это так:
always @(posedge clk) begin
if (!nRst) begin
cordic_angle <= 17'd0;
cordic_quadrant <= 2'd0;
rotator_state <= 2'd0;
end else begin
if (frame_changed) begin
case (rotator_state)
2'd0: begin
if (cordic_angle[15:8] == 8'd89) begin
cordic_quadrant <= cordic_quadrant + 1'b1;
rotator_state <= 2'd1;
end else
cordic_angle[15:8] <= cordic_angle[15:8] + 1'b1;
end
2'd1: begin
if (cordic_angle[15:8] == 8'd1) begin
cordic_quadrant <= cordic_quadrant + 1'b1;
rotator_state <= 2'd0;
end else
cordic_angle[15:8] <= cordic_angle[15:8] - 1'b1;
end
default: rotator_state <= 2'd0;
endcase
end
end
end
Далее значение угла подаётся на модуль CORDIC, который и вычисляет нам значения sin и cos.
cordic CORDIC(
.clk(clk),
.rst(~nRst),
.x_i(17'd19896),
.y_i(16'd0),
.theta_i(cordic_angle),
.x_o(COS),
.y_o(SIN),
.theta_o(),
.valid_in(),
.valid_out()
);
Далее не сложно догадаться, что расчёт координат каждого последующего пикселя будет производиться по формуле:
x’ = cos(angle) * x — sin(angle) * y;
y’ = sin(angle) * x + cos(angle) * y;
Если оставить всё в таком виде, то вращение будет с центром в начале координат. Такое вращение нас не устраивает, нам нужно чтобы картинка вращалась вокруг своей оси с центром в середине изображения. Для этого нам надо вести вычисления относительно центра изображения.
parameter PRECISION = 15;
parameter OUTPUT = 12;
parameter INPUT = 12;
parameter OUT_SIZE = PRECISION + OUTPUT;
parameter BUS_MSB = OUT_SIZE + 2;
wire [15:0] res_x = RES_X - 1'b1;
wire [15:0] res_y = RES_Y - 1'b1;
assign dx = {1'b0, RES_X[11:1]};
assign dy = {1'b0, RES_Y[11:1]};
always @(posedge clk) begin
delta_x <= dx << PRECISION;
delta_y <= dy << PRECISION;
еnd
Далее вычисляем значения cos(angle) * x, sin(angle) * x, cos(angle) * y, sin(angle) * y.
Можно вычислять и так:
always @(posedge clk) begin
mult_xcos <= (xi - dx) * COS;
mult_xsin <= (xi - dx) * SIN;
mult_ycos <= (yi - dy) * COS;
mult_ysin <= (yi - dy) * SIN;
end
Но я решил использовать мегафункции lpm_mult. Их использование значительно повышает Fmax.
reg signed [BUS_MSB: 0] tmp_x, tmp_y, mult_xsin, mult_xcos, mult_ysin, mult_ycos;
reg signed [BUS_MSB: 0] delta_x = 0, delta_y = 0;
wire signed [11:0] dx, dy;
reg signed [BUS_MSB: 0] mxsin, mxcos, mysin, mycos;
reg signed [11:0] ddx, ddy;
always @(posedge clk) begin
ddx <= xi - dx;
ddy <= yi - dy;
end
wire signed [BUS_MSB-1: 0] mult_xcos1;
wire signed [BUS_MSB-1: 0] mult_xsin1;
wire signed [BUS_MSB-1: 0] mult_ycos1;
wire signed [BUS_MSB-1: 0] mult_ysin1;
lpm_mult M1(.clock(clk), .dataa(COS), .datab(ddx), .result(mult_xcos1));
defparam M1.lpm_widtha = 17;
defparam M1.lpm_widthb = 12;
defparam M1.lpm_pipeline = 1;
defparam M1.lpm_representation = "SIGNED";
lpm_mult M2(.clock(clk), .dataa(SIN), .datab(ddx), .result(mult_xsin1));
defparam M2.lpm_widtha = 17;
defparam M2.lpm_widthb = 12;
defparam M2.lpm_pipeline = 1;
defparam M2.lpm_representation = "SIGNED";
lpm_mult M3(.clock(clk), .dataa(COS), .datab(ddy), .result(mult_ycos1));
defparam M3.lpm_widtha = 17;
defparam M3.lpm_widthb = 12;
defparam M3.lpm_pipeline = 1;
defparam M3.lpm_representation = "SIGNED";
lpm_mult M4(.clock(clk), .dataa(SIN), .datab(ddy), .result(mult_ysin1));
defparam M4.lpm_widtha = 17;
defparam M4.lpm_widthb = 12;
defparam M4.lpm_pipeline = 1;
defparam M4.lpm_representation = "SIGNED";
После умножения получаем произведения, знак которых нам необходимо менять в каждом следующем квадранте:
always @(posedge clk) begin
mxcos <= mult_xcos1;
mxsin <= mult_xsin1;
mycos <= mult_ycos1;
mysin <= mult_ysin1;
case (cordic_quadrant)
2'd0: begin
mxsin <= -mult_xsin1;
end
2'd1: begin
mxcos <= -mult_xcos1;
mxsin <= -mult_xsin1;
mycos <= -mult_ycos1;
end
2'd2: begin
mxcos <= -mult_xcos1;
mysin <= -mult_ysin1;
mycos <= -mult_ycos1;
end
2'd3: begin
mysin <= -mult_ysin1;
end
endcase
end
Теперь дело осталось за малым — вычислить сами координаты пикселя:
/*
I II III IV
+ + + - - - - -
+ - + + + - - +
*/
always @(posedge clk) begin
tmp_x <= delta_x + mxcos + mysin;
tmp_y <= delta_y + mycos + mxsin;
end
wire [15:0] xo = tmp_x[BUS_MSB] ? 12'd0: tmp_x[OUT_SIZE-1:PRECISION];
wire [15:0] yo = tmp_y[BUS_MSB] ? 12'd0: tmp_y[OUT_SIZE-1:PRECISION];
Отсекаем пиксели, выходящие за границы изображения:
wire [11:0] xo_t = (xo[11:0] > res_x[11:0]) ? 12'd0 : xo[11:0];
wire [11:0] yo_t = (yo[11:0] > res_y[11:0]) ? 12'd0 : yo[11:0];
И его адрес в памяти:
//addr_out <= yo[11:0] * RES_X + xo[11:0];
И снова используем lpm_mult:
reg [11:0] xo_r, yo_r;
always @(posedge clk) begin
xo_r <= xo_t;
yo_r <= yo_t;
end
wire [28:0] result;
lpm_mult M5(.clock(clk), .dataa(RES_X[11:0]), .datab(yo_r[11:0]), .result(result));
defparam M5.lpm_widtha = 12;
defparam M5.lpm_widthb = 12;
defparam M5.lpm_pipeline = 1;
defparam M5.lpm_representation = "UNSIGNED";
always @(posedge clk) addr_out <= result[22:0] + xo_r[11:0];
Вот, собственно, и всё!
Проблемы метода
Как я уже упоминал выше, данный подход имеет много недостатков. Из-за погрешности вычисления в выходной картинке появляются дыры, чем больше угол поворота, тем больше дыр. Это ещё происходит и по тому, что размеры новой картинки больше чем у оригинала. Этот эффект завётся aliasing и существуют методы борьбы с ним, например, медианный фильтр, расмотренный в моей предыдущей статье.
Перед каждым последующим кадром не мешало бы почистить память от предыдущего кадра, чтобы новое изображение получалось на чистом фоне, но это занимает время и придётся пропускать один кадр.
Единственным достоинством метода является простота реализации и скорость обработки т.к. координаты вычисляются на-лету.
Вот что из этого получилось
Ссылки по теме
> CORDIC на русском
> CORDIC for dummies
> CORDIC FAQ
Архив проекта в Квартусе
> Ссылка на яндекс диск.
Комментарии (32)
FGV
30.03.2017 20:48+2Можно же не считать синусы косинусы в зависимости от угла «на лету», достаточно задать начальные значения для нулевого угла, и две константы — sin(1градус), cos(1градус).
Последующие значения вычисляются по формулам суммы углов синусов и косинусов, в итоге для расчета синуса/косинуса понадобится два умножения и одно сложение.
leshabirukov
31.03.2017 14:59Даже формула суммы углов не нужна, сам CORDIC позволяет легко доворачивать на малый угол, если использовать результат с предыдущего шага.
Algoritmist
31.03.2017 15:03+1Здесь синус и косинус угла надо вычислять заранее только один раз на кадр. Как вариант можно таблицу записать и интерполировать. Но зачем это делать, если это уже сделано в CORDIC?
При рекуррентном вычислении надо быть осторожным. Ошибка копиться будет. Можно обратно в 0 градусов не вернуться.leshabirukov
31.03.2017 21:28Вот как раз проходя мимо нуля ошибку и можно сбрасывать — в этот момент считаем не от предыдущего, а от нуля.
FGV
01.04.2017 21:25Но зачем это делать, если это уже сделано в CORDIC?
Да он жрет много ресурсов, и если без него можно обойтись, нафиг он нужен?
При рекуррентном вычислении надо быть осторожным. Ошибка копиться будет. Можно обратно в 0 градусов не вернуться.
Более того, не вернется в 0 гарантировано, ошибка будет только расти. Выход один — обнуление при проходе значений 90*n градусов.
Foveator
30.03.2017 22:29Мне было бы интересным реализация вращения видеоизображения без многократных вычислений.
Например коллега Вячеслав подписался, что в ПЛИС реализовал алгоритм VGA
https://habrahabr.ru/post/157863/
А где VGA там и RGB разъёмы. А там и сигналы строчной и кадровой развертки идут по двум отдельным каналам.
Предлагаю поворачивать видеоизображения управлением кадровой и строчной развертки.
Готов участвовать в написании алгоритма (правда в Фортране).
ubobrov
30.03.2017 22:35У меня HDMI.
Боюсь оказаться неправым, но что-то мне подсказывает, что это не совсем возможно.Foveator
31.03.2017 11:36Душа боится, а руки делают!
Вам прислать переходник HDMI-VGA. Может у Вас тамперированный монитор поэтому в нём нет CRT 640*480 и он не позволяет по сигналам строчной и кадровой развертки dx и dy выводить пиксель в заданную точку экрана?
Полагаю, что повесив в подъезде объявление «Дорогие соседи! Кому не жалко, занесите свой монитор VGA на пятый этаж он мне нужен.» решит проблему.
При реализации простого алгоритма надо следить, чтобы края вращаемого видеоизображения в мониторе не вылезли за края видеоизображения, формируемого видеокамерой. Поэтому видеокадр в видеокамере должен иметь размерность не менее 800 пикселей по каждой стороне.
Я весь внимания и в готовности к сотрудничеству советами и средствами.ubobrov
31.03.2017 11:46Не, переходника мне не надо. У меня выход HDMI и вход на мониторе HDMI, а VGA ну совсем нет ни где.
Foveator
31.03.2017 13:19Не хотите со стороны монитора делать моментальное преобразование видеоизображение.
Попробуем подойти к проблеме со стороны целесообразности работы.
Предположим, что камера у Вас на подвижном основании снимает улицу и при этом качается.
С камерой качается гироскоп и выдает оперативно информацию о наклоне камеры.
Сигнал поступает к Вам на FPGA.
Камера и гироскоп можно использовать, установленные в смартфоне.
Вы, используя свой алгоритм, поворачиваете изображение в сторону. противоположную наклону камеры, чтобы для зрителя изображение улицы не качалось.
Задача компенсации вибрации камеры кажется может быть востребована.ubobrov
31.03.2017 13:26Именно эта идея мне и приходила в голову ранее. Да, задача интересная и я ей обязательно займусь, но на очереди у меня blob detection с последующим распознованием на нейронных сетях.
oleg1977
30.03.2017 23:30Не знаю реализацию, которой воспользовались вы, но оригинальный кордик сразу на выходе дает повернуты вектор, умножения и сложения делать не надо. Вычисление sin и cos производятся поворотом вектора (1.0; 0.0), если вместо него подать на вход кордика (x;y) то на выходе сразу будет повернутый вектор
antonsosnitzkij
31.03.2017 00:02+1Не знаю что Вы подразумеваете под оригинальным кордиком, но настоящий кордик поворачивает вектор и далее в зависимости от знаков делает сложение умножение/сдвиг. Сам по себе поворот осуществляется с помощью сложения и умножения
oleg1977
31.03.2017 15:10Классический кордик описан, например, здесь
Он преобразует (x0; y0; z0) --> K(x0*cos(z0)-y0*sin(z0); x0*sin(z0)+y0*cos(z0); 0)
Отсюда, если взять x0 = 1/K, y0=0, то получим (1/K; 0; z0) --> (cos(z0); sin(z0); 0)
То есть, вычисление косинусов и синусов — частный случай оригинального кордика.
Я лишь обращал внимание, что вместо преобразований:
1) вычислить cos(z0), вычислить sin(z0)
2) вычислить x0*cos(z0)-y0*sin(z0)
3) вычислить x0*sin(z0)+y0*cos(z0)
можно обойтись одним:
1) вычислить x0*cos(z0)-y0*sin(z0) и x0*sin(z0)+y0*cos(z0)oleg1977
31.03.2017 15:16Далее значение угла подаётся на модуль CORDIC, который и вычисляет нам значения sin и cos.
cordic CORDIC( .clk(clk), .rst(~nRst), .x_i(17'd19896), .y_i(16'd0), .theta_i(cordic_angle), .x_o(COS), .y_o(SIN), .theta_o(), .valid_in(), .valid_out() );
Далее не сложно догадаться, что расчёт координат каждого последующего пикселя будет производиться по формуле:
x’ = cos(angle) * x — sin(angle) * y;
y’ = sin(angle) * x + cos(angle) * y;
Вместо этого можно было бы написать что-то вроде:
cordic CORDIC( .clk(clk), .rst(~nRst), .x_i(x), .y_i(y), .theta_i(cordic_angle), .x_o(x'), .y_o(y'), .theta_o(), .valid_in(), .valid_out() );
vassabi
30.03.2017 23:37второй вариант — это использовать честную координату (ну, почти честную) для расчета начала строки, а для следующих — алгоритм брезенхема.
Брезенхем же в чем — есть виртуальная шестеренка с К зубцами (если использовать вашу КПДВ и рисовать строки слева направо от 0,0 до W,0 и снизу вверх), на каждый ее шаг «поворота» увеличиваем координату х на 1, а как сделает полный оборот — увеличиваем координату у на 1, т.д — по крайней мере до угла 45 градусов :)
Если больше 45 — то там по столбцам нужно идти. Но идея от єтого не меняется.
При угле поворота картинки 0 градусов — К равно W (т.е. ни разу у не увеличивается), при 45 — К равно 1 (т.е. у увелчивается так же как и х — а каждом шаге).
Так как строчки идут параллельно, то у всех у них будет одинаковый угол — и параметр К. Отличаться только будет координата начала и фаза виртуальной «шестерни». Так что нужно будет честно вычислить только в начале кадра — куда придет 0,0 и какая в этой точке фаза (интуитивно кажется — что можно задать равной 0 — типа начало, но хз, надо нарисовать или попробовать).
После этого — для каждой следуюущей строки можно вычислять только стартовую х_і координату (а у — получить из брезенхема, так как это y_0 + i — (х_0 — х_і) / К. Заодно и дырок между строчками не будет). Начальная фаза при этом будет равна (х_0 — х_і) mod К.
Далее — можно заметить еще, что если размеры (H, W) исходного и результирующего прямоугольников равны, и мы делаем поворот вокруг цента, то если мы сейчас читаем из точки с координатой (х_і, у_і), и после преобразования, ее координаты (х2_і, у2_і) находятся вне результирующего прямоугольника, то это значит, что в результирующем прямоугольнике есть точка (х_і, H — у_і), у которой нужно выставить цвет фона :)
(В принципе, это соображение будет работать при любом угле поворота… вроде бы ...)
если это все можно выразить на верилоге — было бы интересно посмотреть результат :DDDubobrov
30.03.2017 23:59Не уверен, что у меня получится это реализовать )).
Там есть ещё одна проблема, которая заключается в скорости чтения/записи данных во внешнюю память. Сейчас камера выдаёт 25 fps. С разрешением 320x240 времени хватает на запись данных не бёрстами. При увеличении угла поворота возможность записи бёрстами сходит на нет т.к. адреса пикселей находятся не в одном ряду, а это означает что каждую запись надо делать начиная с открытия ряда, выдержки tRRD, tRP для SDRAM, закрытия ряда. Это очень сильно влияет на производительность. Отсюда, если выбирать данные из SDRAM по столбцам, то на каждый пиксель надо тратить уйму времени, да плюс ещё и входные данные от камеры надо успевать записывать, а HDMI контроллер не ждёт, он работает на 25Мгц. По этому для данной реализации куда дешевле по времени все процедуры преобразования производить на записи, а чтение пусть работает бёрстами по строкам.vassabi
31.03.2017 02:22хмм, а возвращаясь к идее из первого комментария:
если в памяти буферы расположить не подряд — А0A1...Am, B0B1...Bm, C0C1...Cm, а чередуя, чтобы А0B0C0[ВЫРАВНИВАНИЕ], А1B1C1[ВЫРАВНИВАНИЕ]… АmBmCm[ВЫРАВНИВАНИЕ],
и потом на приходе очередного сигнала вычислять координату чтения, читать данные из промежуточного буфера, а потом делать запись данных из камеры на рабочий и запись прочитанного (или черного/фонового) цвета на выходной буфер делать сразу одним заходом (только надо, чтобы коррекция адресов для буферов А/B/C происходила прямо в проводах)
так можно же?
Algoritmist
31.03.2017 11:11+2Поворот на ПЛИС можно сделать двумя последовательными косыми сдвигами. Отдельно по строкам, отдельно по столбцам. Тогда с памятью идет поточная работа. В результате выходит последовательное применение двух одномерных фильтров. Справедливости ради, стоит сказать, что это не совсем «косой сдвиг», но суть близкая. Пусть сначала идет обработка по строкам (матрица A), а потом по столбцам (матрица B). Матрица поворота (М) раскладывается так (для аффинных преобразований такой прием проходит):
B * A = M
A вида (косой горизонтальный сдвиг)
a1 a2
1 0
B вида (косой вертикальный сдвиг)
1 0
b1 b2
т.е.
x* = a1 x + a2 y
y* = b1 x* + b2 y (стоит обратить внимание, что «x» тут новое)
Осталось найти коэффициенты (a1, a2, b1, b2) этих матриц. Выходит система уравнений. Ответ:
a1 = cos(alpha)
a2 = sin(alpha)
b1 = cos(alpha) / sin(alpha) = ctg(alpha)
b2 = -1 / sin(alpha)Algoritmist
31.03.2017 15:33Извините, в матрице A ошибка:
a1 a2
0 1
и ответ:
a1 = cos(alpha)
a2 = sin(alpha)
b1 = -sin(alpha) / cos(alpha) = -tg(alpha)
b2 = 1 / cos(alpha)
Еще такой подход не работает на углах поворота близких к 90 градусам. Можно работать от -45 до 45, а иначе работать с модификацией, которая как бы поворачивающей заранее изображение на 90 градусов (меняет х и у местами).
msts2017
31.03.2017 14:47а на этом железе есть какая-то разница целочисленные вычисления или нет?
Algoritmist
31.03.2017 14:55Есть. Для простоты можно считать, что числа только целочисленные (с фиксированной точностью).
vassabi
ничего не знаю в верилоге — а вы можете не записывать, а вычитывать данные по заданному адресу?
Тогда можно было бы рисовать тем же формулами, но с вращением координат чтения — что решало бы и проблемы с пропусками точек и с фоном, когда нет данных по краям :)
Ну или сделать конвейер — выделить два массива (рабочий=0 и буфер=1) и текущую точку просто писать в буфер, а на экран в эту же точку выводить из повернутых координат буфера (который был прошлым кадром). После прохода — поменять ссылки на них (рабочий=1 а буфер=0) %)
ubobrov
Ну это на Си так легко сделать, на ПЛИС не знаю. Вы предлагаете сначала записать кадр без изменения, затем расчитать параметры нового изображения, затем вычитывать из записанного кадра в новый расчитанный кадр, а потом снова вычитывать из нового и выводить на экран. Этот метод хорош при наличии фреймбуфера, который заменит нам «буфер=1», у меня его нет.
Там проще можно сделать. Можно ресайз нового кадра делать чтобы он вписывался в наш фрейм 320x240.
Дело ещё и в том, что уменя места в ПЛИС на подобную обвязку может нехватить. Сам CORDIC отожрал неприлично много ресурсов.
bambruysk
На ПЛИС это делается легко. Для каждого буфера добавляется базовый адрес. Его изменение и дает переключение буферов. В целом, на каком железе это работает? Фреймбуфер во внешней памяти или в блочной памяти плисины?
ubobrov
А нет никакого фреймбуфера, есть только область памяти в SDRAM в которую пишется кадр, а затем читается из неё в FIFO, а оттуда в HDMI. Работает на EP2C8F256.
Идею я понял. Попробую повращать читателя, а не писателя как сейчас сделано.
bambruysk
Я это имел в виду. А вы пробовали вычислять синусы и косинусы табличным методом?
ubobrov
Нет. Мне интересен был именно CORDIC. На него у меня есть планы в дальнейших разработках.