{
作者 ： 刘留
参考文献为 ：
Жан-Ив Бугит "Пирамидальная реализация Лукаса Канаде Описание Алгоритма" Лукас Канаде Описание алгоритма "
http://www.aivisoft.net/
Geo.cra [at] gmail [dot] com
}
единица OpticalFlowlk;
интерфейс
Использование
Математика, окна, Sysutils, варианты, классы, графика, Unitypes, ColorConv;
тип
TopticalFlowlk = Class
Частный
Imageold, ImageNew: ttripleLongIntArray;
ImageGray, DX, DY, DXY: TdoubleLongintArray;
Собственные значения: tdoubleextendendarray;
Wbpoint: tdoublebooleanarray;
Высота, ширина, L, Radiox, Radioy: Longint;
Процедура CornerDetect (Swidth, Sheuett: Longint; качество: расширенное);
Процедура MakePyramid (var Images: ttriplelongintarray; Swidth, Sheuett, SL: Longint);
публичный
Кадр: tbitmap;
Особенности: tsinglepointinfoarray;
FeatureCount, SupportCount: Longint;
Скорость, радио: расширенная;
Процедура init (Swidth, Sheuett, SL: Longint);
Процедура initfeatures (кадры: tbitmap);
Процедура CalopticalFlowlk (кадры: tbitmap);
разрушитель разрушил; переопределить;
конец;
выполнение
Процедура TopticalFlowlk.cornerDetect (Swidth, Sheuet: Longint; качество: расширенное);
вар
Я, J, Fi, FJ: Longint;
A, B, C, Sum, Minaccept, MaxeigenValue: Extended;
начинать
FeatureCount: = 0;
{
下面采用 Хорошая функция для отслеживания 介绍的方法
J. Shi и C. Tomasi «Хорошие функции для отслеживания», CVPR 94
}
для i: = 1 к Swidth - 2 do
для j: = 1 к шейту - 2 начинать
dx [i, j]: = imagegray [i - 1, j - 1] + 2 * imagegray [i - 1, j] + imagegray [i - 1, j + 1]
- (Imagegray [i + 1, j - 1] + 2 * Imagegray [i + 1, j] + imagegray [i + 1, j + 1]);
dy [i, j]: = imagegray [i - 1, j + 1] + 2 * imagegray [i, j + 1] + imagegray [i + 1, j + 1]
- (Imagegray [i - 1, j - 1] + 2 * imagegray [i, j - 1] + imagegray [i + 1, j - 1]);
dxy [i, j]: = imagegray [i + 1, j - 1] + imagegray [i - 1, j + 1]
- (Imagegray [i - 1, j - 1] + imagegray [i + 1, j + 1]);
конец;
{求取 Sobel 算子的 dx, dy, dxy
DX:
| 1 0 -1 |
| 2 0 -2 |
| 1 0 -1 |
Dy:
| -1 -2 -1 |
| 0 0 0 |
| 1 2 1 |
Dxy
| -1 0 1 |
| 0 0 0 |
| 1 0 -1 |}
Maxeigenvalue: = 0;
для i: = 2 к Swidth - 3
для j: = 2 к шейту - 3 начинать
A: = 0; b: = 0; c: = 0;
для fi: = i - 1 до i + 1 do
для fj: = j - 1 до j + 1 do begin
a: = a + sqr (dx [fi, fj]);
b: = b + dxy [fi, fj];
c: = c + sqr (dy [fi, fj]);
конец;
A: = A / 2; C: = C / 2;
Собственные значения [i, j]: = (a + c - sqrt ((a - c) * (a - c) + b * b));
Если собственные значения [i, j]> maxeigenvalue, то maxeigenvalue: = собственные значения [i, j];
конец;
{求取矩阵
| ∑dx*dx ∑dxy |
M = | |
| ∑dxy ∑dy*dy |
的特征值
λ = ∑dx*dx + ∑dy*dy - ((∑dx*dx + ∑dy*dy)^2-4*(∑dx*dx*∑dy*dy - ∑dxy*∑dxy))^1/2 }
Minaccept: = maxeigenvalue * качество;
{设置最小允许阀值}
для i: = 8 до Swidth - 9
для j: = 8 к шейту - 9
Если собственные значения [i, j]> minaccep
Wbpoint [i, j]: = true;
Inc (featureCount);
конец еще
Wbpoint [i, j]: = false;
для i: = 8 до Swidth - 9
для j: = 8 к шейту - 9
Если wbpoint [i, j], то начните
sum: = собственные значения [i, j];
для fi: = i - 8 до i + 8, начинайте
для fj: = j - 8 до j + 8 do
Если sqr (fi - i) + sqr (fj - j) <= 64, тогда
if (собственные значения [fi, fj]> = sum) и ((fi <> i) или (fj <> j)) и (wbpoint [fi, fj]), затем начните
Wbpoint [i, j]: = false;
Дек (featureCount);
перерыв;
конец;
Если не wbpoint [i, j], то сломайте;
конец;
конец;
{用非最大化抑制来抑制假角点}
setLength (функции, featureCount); fi: = 0;
для i: = 8 до Swidth - 9
для j: = 8 к шейту - 9
Если wbpoint [i, j], то начните
Особенности [fi] .info.x: = i;
Особенности [fi] .info.y: = j;
Особенности [fi] .index: = 0;
Inc (Fi);
конец;
{输出最终的点序列}
конец;
Процедура TopticalFlowlk.init (Swidth, Sheuth, SL: Longint);
начинать
Ширина: = Swidth; Высота: = Шейт; L: = sl;
SetLength (ImageOLD, ширина, высота, L);
setLength (ImageNew, ширина, высота, L);
Кадр: = tbitmap.create;
Frame.width: = ширина; Frame.height: = высота;
Frame.pixelformat: = pf24bit;
SetLength (ImageGray, ширина, высота);
SetLength (собственные значения, ширина, высота);
setLength (dx, ширина, высота);
setLength (dy, ширина, высота);
setLength (dxy, ширина, высота);
SetLength (WBPoint, ширина, высота);
FeatureCount: = 0;
конец;
Процедура TopticalFlowlk.makePyramid (var Images: ttriplelongintarray; Swidth, Sheuett, SL: Longint);
вар
I, J, K, II, JJ, Nwidth, Nheight, Owidth, Oheight: Longint;
начинать
{生成金字塔图像}
Owidth: = Swidth; oheight: = sheueth;
для k: = 1 до Sl - 1
nwidth: = (Owidth + 1) shr 1; nheight: = (oheight + 1) shr 1;
Для i: = 1 к nwidth - 2 начинать
ii: = i shl 1;
для j: = 1 до nheight - 2 начинать
JJ: = J Shl 1;
Изображения [i, j, k]: = (Image [ii, jj, k - 1] shl 2 +
Изображения [II - 1, JJ, K - 1] SHL 1 + изображения [II + 1, JJ, K - 1] SHL 1 + изображения [II, JJ - 1, K - 1] SHL 1 + Изображения [II, JJ + 1, k - 1] shl 1 +
Изображения [ii - 1, jj - 1, k - 1] + изображения [ii + 1, jj - 1, k - 1] + изображения [II - 1, JJ + 1, k - 1] + изображения [ii + 1 , jj + 1, k - 1]) shr 4;
{高斯原则 ， shl 右移位 ， shr 左移位}
конец;
конец;
Для i: = 1 к nwidth - 2 начинать
ii: = i shl 1;
Изображения [i, 0, k]: = (изображения [II, 0, K - 1] Shl 2 +
Изображения [ii - 1, 0, k - 1] Shl 1 + изображения [ii + 1, 0, k - 1] Shl 1 + изображения [II, 0, K - 1] Shl 1 + изображения [II, 1, k - 1] SHL 1 +
Изображения [ii - 1, 0, k - 1] + изображения [ii + 1, 0, k - 1] + изображения [ii - 1, 1, k - 1] + изображения [ii + 1, 1, k - 1 ]) SHR 4;
Изображения [i, nheight - 1, k]: = (Images [ii, oheight - 1, k - 1] shl 2 +
Изображения [ii - 1, oheight - 1, k - 1] shl 1 + изображения [ii + 1, oheight - 1, k - 1] shl 1 + изображения [ii, oheight - 2, k - 1] shl 1 + изображения [ii, oheight - 1, k - 1] shl 1 +
Изображения [ii - 1, oheight - 2, k - 1] + изображения [ii + 1, oheight - 2, k - 1] + изображения [ii - 1, oheight - 1, k - 1] + изображения [ii + 1 , oheight - 1, k - 1]) shr 4;
конец;
{处理上下边}
для j: = 1 до nheight - 2 начинать
JJ: = J Shl 1;
Изображения [0, j, k]: = (изображения [0, JJ, K - 1] Shl 2 +
Изображения [0, JJ, K - 1] Shl 1 + изображения [1, JJ, K - 1] Shl 1 + изображения [0, JJ - 1, K - 1] Shl 1 + изображения [0, JJ + 1, k - 1] SHL 1 +
Изображения [0, JJ - 1, K - 1] + изображения [1, JJ - 1, K - 1] + изображения [0, JJ + 1, K - 1] + изображения [1, JJ + 1, k - 1 ]) SHR 4;
Изображения [nwidth - 1, j, k]: = (изображения [Owidth - 1, JJ, K - 1] Shl 2 +
Изображения [Owidth - 2, JJ, K - 1] Shl 1 + изображения [Owidth - 1, JJ, K - 1] Shl 1 + изображения [Owidth - 1, JJ - 1, K - 1] Shl 1 + Изображения [Owidth - 1, jj + 1, k - 1] shl 1 +
Изображения [Owidth - 2, JJ - 1, K - 1] + изображения [Owidth - 1, JJ - 1, K - 1] + изображения [Owidth - 2, JJ + 1, k - 1] + изображения [Owidth - 1 , jj + 1, k - 1]) shr 4;
конец;
{处理左右边}
Изображения [0, 0, k]: = (изображения [0, 0, k - 1] Shl 2 +
Изображения [0, 0, k - 1] Shl 1 + изображения [1, 0, K - 1] Shl 1 + изображения [0, 0, K - 1] Shl 1 + изображения [0, 1, k - 1] Shl 1 +
Изображения [0, 0, K - 1] + изображения [1, 0, K - 1] + изображения [0, 1, K - 1] + изображения [1, 1, K - 1]) SHR 4;
{处理左上点}
Изображения [0, nheight - 1, k]: = (изображения [0, oheight - 1, k - 1] shl 2 +
Изображения [0, oheight - 1, k - 1] shl 1 + изображения [1, oheight - 1, k - 1] shl 1 + изображения [0, oheight - 2, k - 1] shl 1 + изображения [0, Oheight - 1, k - 1] Shl 1 +
Изображения [0, oheight - 2, k - 1] + изображения [1, oheight - 2, k - 1] + изображения [0, oheight - 1, k - 1] + изображения [1, oheight - 1, k - 1 ]) SHR 4;
{处理左下点}
Изображения [nwidth - 1, 0, k]: = (изображения [Owidth - 1, 0, k - 1] Shl 2 +
Изображения [Owidth - 2, oheight - 1, k - 1] shl 1 + изображения [Owidth - 1, oheight - 1, k - 1] Shl 1 + изображения [Owidth - 1, oheight - 1, k - 1] Shl 1 + Изображения [Owidth - 1, oheight - 1, k - 1] shl 1 +
Изображения [Owidth - 2, oheight - 1, k - 1] + изображения [Owidth - 1, oheight - 1, k - 1] + изображения [Owidth - 2, oheight - 1, k - 1] + изображения [Owidth - 1 , oheight - 1, k - 1]) shr 4;
{处理右上点}
Изображения [nwidth - 1, nheight - 1, k]: = (изображения [Owidth - 1, oheight - 1, k - 1] Shl 2 +
Изображения [Owidth - 2, oheight - 1, k - 1] Shl 1 + изображения [Owidth - 1, Oheight - 1, K - 1] Shl 1 + изображения [Owidth - 1, Oheight - 2, k - 1] Shl 1 + Изображения [Owidth - 1, oheight - 1, k - 1] shl 1 +
Изображения [Owidth - 2, oheight - 2, k - 1] + изображения [Owidth - 1, oheight - 2, k - 1] + изображения [Owidth - 2, oheight - 1, k - 1] + изображения [Owidth - 1 , oheight - 1, k - 1]) shr 4;
{处理右下点}
конец;
конец;
Процедура TopticalFlowlk.InitFeatures (рамки: TBITMAP);
вар
Я, J: Longint;
Линия: prgbtriple;
начинать
SetStrechBltMode (frame.canvas.handle, strept_halftone);
Strectblt (frame.canvas.handle, 0, 0, ширина, высота, кадры. Canvas.handle, 0, 0, Frames.width, Frames.height, Srccopy);
для i: = 0 до высоты - 1 дела
Line: = Frame.Scanline [i];
для j: = 0 до ширины - 1 дела
ImageOLD [j, i, 0]: = (Line^.rgbtblue * 11 + line^.rgbtgreen * 59 + line^.rgbtred * 30) div 100;
Imagegray [j, i]: = Imageold [j, i, 0];
Inc (Line);
конец;
конец;
{初始化金字塔图像第一层 Imageold [x, y, 0]}
MakePyramid (Imageold, ширина, высота, L);
{生成金字塔图像}
CornerDetect (ширина, высота, 0,01);
{进行强角点检测}
конец;
процедура topticalFlowlk.calopticalflowlowk (рамки: tbitmap);
вар
I, J, Fi, FJ, K, LL, M, DX, DY, GX, GY, PX, PY, KX, KY, ED, EDC, NWIDTH, Nheight: Longint;
NX, NY, VX, VY, A, B, C, D, E, F, IK: Extended;
Ix, iy: tdoubleextendendarray;
Линия: prgbtriple;
Изменение: логический;
начинать
SetStrechBltMode (frame.canvas.handle, strept_halftone);
Strectblt (frame.canvas.handle, 0, 0, ширина, высота, кадры. Canvas.handle, 0, 0, Frames.width, Frames.height, Srccopy);
для i: = 0 до высоты - 1 дела
Line: = Frame.Scanline [i];
для j: = 0 до ширины - 1 дела
ImageNew [j, i, 0]: = (line^.rgbtblue * 11 + line^.rgbtgreen * 59 + line^.rgbtred * 30) div 100;
Inc (Line);
конец;
конец;
{初始化金字塔图像第一层 ImageNew [x, y, 0]}
MakePyramid (ImageNew, ширина, высота, L);
{生成金字塔图像}
SetLength (IX, 15, 15); SetLength (IY, 15, 15);
{申请差分图像临时数组}
для m: = 0 к featurecount - 1 do begin
{算法细节见:
Жан-Ив Бугет "Пирамидальная реализация Лукаса Канаде Описание алгоритма"}}}
gx: = 0; GY: = 0;
для ll: = l - 1 до 0
px: = функции [m] .info.x shr ll;
py: = особенности [m] .info.y shr ll;
{对应当前金字塔图像的 u 点 ： u [l]: = u/2^l}
nwidth: = ширина shr ll; nheight: = высота shr ll;
A: = 0; B: = 0; C: = 0;
для i: = max (px - 7, 1) до мин (px + 7, nwidth - 2) do
для j: = max (py - 7, 1) до минима
fi: = i - px + 7; fj: = j - py + 7;
Ix [fi, fj]: = (Imageold [i + 1, j, ll] - Imageold [i - 1, j, ll]) / 2;
Iy [fi, fj]: = (Imageold [i, j + 1, ll] - Imageold [i, j - 1, ll]) / 2;
A: = a + ix [fi, fj] * ix [fi, fj]; B: = b + ix [fi, fj] * iy [fi, fj];
C: = c + iy [fi, fj] * iy [fi, fj];
конец;
{计算 2 阶矩阵 g:
| Ix (x, y)*ix (x, y) ix (x, y)*iy (x, y) |
∑ | ix (x, y)*iy (x, y) iy (x, y)*iy (x, y) |}
D: = a * c - b * b;
vx: = 0; vy: = 0; dx: = 0; dy: = 0;
Если ABS (D)> 1E-8, начните
для k: = 1-10
E: = 0; F: = 0;
для i: = max (px - 7, 1) до мин (px + 7, nwidth - 2) do
для j: = max (py - 7, 1) до минима
fi: = i - px + 7; fj: = j - py + 7;
kx: = i + gx + dx; ky: = j + gy + dy;
Если kx <0, то kx: = 0; Если kx> nwidth - 1, то kx: = nwidth - 1;
Если ky <0, то ky: = 0; Если ky> nheight - 1, то ky: = nheight - 1;
IK: = ImageOLD [i, J, Ll] - ImageNew [KX, KY, LL];
E: = e + ik * ix [fi, fj];
F: = f + ik * iy [fi, fj];
конец;
{计算 2x1 阶矩阵 b:
| Ik (x, y)*ix (x, y) |
∑ | ik (x, y)*iy (x, y) |}
nx: = (c * e - b * f) / d;
NY: = (b * e - a * f) / (-d);
{计算 η = g^-1*b}
vx: = vx + nx; vy: = vy + ny;
dx: = trunc (vx); dy: = tunc (vy);
{得到相对运动向量 d}
конец;
конец;
gx: = (gx + dx) shl 1; GY: = (GY + DY) SHL 1;
{得到下一层的预测运动向量 g}
конец;
GX: = GX Div 2; GY: = Gy Div 2;
px: = px + gx; py: = py + gy;
Особенности [m] .info.x: = px;
Особенности [m] .info.y: = py;
Особенности [m] .vector.x: = gx;
Особенности [M] .Vector.Y: = GY;
if (px> ширина - 1) или (px <0) или (py> высота - 1) или (py <0), затем функции [m] .index: = 1;
{失去特征点处理}
конец;
для k: = 0 до l - 1, начинается
nwidth: = ширина shr k; nheight: = высота shr k;
для i: = 0 к nwidth - 1 do
для j: = 0 до nheight - 1 do
ImageOLD [i, J, K]: = ImageNew [i, J, K];
конец;
{复制 j 到 i}
повторить
Изменение: = ложь;
для i: = 0 к featureCount - 1 Do Begin
Если функции [i] .index = 1, тогда
для j: = i + 1 для featurecount - 1 do begin
Особенности [j - 1]: = функции [j];
Изменение: = true;
конец;
Если изменить, то перерыв;
конец;
Если изменить, то dec (featureCount);
пока не изменится;
setLength (функции, featureCount);
{删除失去的特征点}
Скорость: = 0; Радио: = 0; Радиосвязь: = 0; Radioy: = 0;
Если featureCount> 0, тогда начинайте
SupportCount: = 0;
для i: = 0 к featureCount - 1 do
if (функции [i] .vector.x <> 0) и (функции [i] .vector.y <> 0), затем начните
Radiox: = radiox + функции [i] .vector.x;
Radioy: = Radioy + особенности [i] .vector.y;
Скорость: = скорость + sqrt (sqr (функции [i] .vector.x) + sqr (функции [i] .vector.y));
Inc (поддержка);
конец;
Если поддержка, затем 0, тогда начните
Скорость: = скорость / поддержка; //*0.0352778;
Радио: = arctan2 (radioy, radiox);
конец;
конец;
{计算平均速度和整体方向}
Frame.canvas.pen.style: = pssolid;
Frame.canvas.pen.width: = 1;
Frame.canvas.pen.color: = cllime;
Frame.canvas.brush.style: = bsclear;
для i: = 0 к featureCount - 1 do
Frame.canvas.ellipse (функции [i] .info.x - 4, функции [i] .info.y - 4, функции [i] .info.x + 4, функции [i] .info.y + 4) ;
{用绿圈标识特征点}
Frame.canvas.pen.color: = clyellow;
для i: = 0 к featureCount - 1 Do Begin
Frame.canvas.moveto (функции [i] .info.x - функции [i] .vector.x, функции [i] .info.y - функции [i] .vector.y);
Frame.canvas.lineto (функции [i] .info.x, функции [i] .info.y);
конец;
{用黄色线条表示运动向量}
if (supportCount> 0) и (скорость> 3), затем начните
Frame.canvas.pen.style: = pssolid;
Frame.canvas.pen.width: = 6;
Frame.canvas.pen.color: = clwhite;
Frame.canvas.ellipse (ширина div 2 - 100, высота div 2 - 100, width div 2 + 100, высота div 2 + 100);
Frame.canvas.moveto (ширина Div 2, высота div 2);
Frame.canvas.pen.color: = clblue;
Frame.canvas.lineto (ширина div 2 + trunc (90 * cos (радио)), высота div 2 + trunc (90 * sin (радио)));
Frame.canvas.pen.style: = psclear;
Frame.canvas.brush.style: = bssolid;
Frame.canvas.brush.color: = clred;
Frame.canvas.ellipse (ширина div 2 + trunc (90 * cos (радио)) - 8, высота div 2 + trunc (90 * sin (радио)) - 8, ширина div 2 + trunc (90 * cos (радио) ) + 8, высота div 2 + trunc (90 * sin (радио)) + 8);
конец;
конец;
деструктор TopticalFlowlk.destroy;
начинать
SetLength (ImageOLD, 0);
SetLength (ImageNew, 0);
SetLength (ImageGray, 0);
SetLength (собственные значения, 0);
setLength (dx, 0);
setLength (dy, 0);
setLength (dxy, 0);
SetLength (WBPoint, 0);
унаследован;
конец;
конец.