一、简介
根据对动态图像和时变图像序列的分析，来确定客观物体与观察者之间的相对运动参数，是当今计算机视觉研究领域中的一个热门课题。在医疗，工业，国防等方面都具有非常重要的现实意义。
本文将介绍采用光流模型分析运动图形：根据图像像素强度守恒原理，建立光流约束方程，计算运动参数，最后结合实例计算两帧样本图像之间的水平和垂直位移量，并绘制光流矢量图。
1 重要概念：光流
光流，是空间运动物体在被观测表面上的像素点运动的瞬时速度场，包含了物体与成像传感器系统之间的相对运动的关系。
2 光流约束方程
物体在空间上一般是相对连续运动，故投影到视网膜/传感器平面上的图像应该也是连续变化的。因此可以假设运动图像函数 ƒ(x,y,t) 是关于变量x,y,t的连续函数。
故，我们可以认为，t 时刻和 t+dt 时刻的图像点强度相等：
然后，上式右边在 (x,y,t) 做泰勒级数展开，约去高阶项并同时除以dt，可得：
这就是光流约束方程。
这个式子表示：图像强度对时间变化率等于强度的空间变化率与运动速度的乘积( 写成向量式更直观 )
记：
代入( 2 )中，可得光流约束方程的第二种形式：
3 光流约束方程的适用性分析
3.1 实际应用中，图像可能有强度不连续的点存在。在这种情况下，只要在不连续区域周围图像强度均匀变化，光流约束方程依然成立。
3.2 需确保图像平面上每一点的变化完全由物体与观察者之间的相对运动引起。
4 光流约束方程的求解
光流约束方程的求解，也即计算两个未知速度分量 u 和 v 。
由 (3) 式知，u 和 v 的解被约束在空间的一条直线上，没有另一个约束条件，我们无法唯一地确定 u 和 v 的解。因此必须附加一个新的约束条件。
显然，u 和 v 随着像素点移动而发生的改变是缓慢的，局部区域的变化不大，尤其是在目标做无变形刚体运动时，局部区域速度的空间变化率为0。
令 uavg 和 vavg 分别表示 u 和 v 领域内速度的平均值，根据上述分析，可知：
这个约束方程反映了速度变化率应该满足的必要条件，而光流约束方程反映了速度与灰度变化率之间的相对关系。
结合光流约束方程和约束条件，对于全部的像素点，可以建立如下的极小化方程：
其中 ε 是极小化问题的总误差，λ 是附加约束的拉格朗日乘子。
为使 ε 最小，令 ε 分别对 u 和 v 求导并取导数为 0 ，得：
可导出：
可用迭代法求解上式，从而得出 u 和 v。
在实际求解中，还需要计算图像强度在空间和时间的偏导数 Ex，Ey，Et 。设ƒi, j, k 为图像在 (x,y,t) 点的强度，ƒi, j,k+1为 (x, y, t+dt) 点的强度，那么 Ex，Ey，Et 可以按照以下式子近似计算出来：
对于全部的 x 和 y 值，当迭代运算的结果满足事先给定的估计容差 ε1 和 ε2 时，迭代计算结束，也即：
5 实际分析步骤
5.1初始化：迭代次数 n = 1，λ = 1，u(n-1) = v(n-1) = 0 。
5.2 根据公式 (9) 计算图像强度的偏导数 Ex，Ey，Et 。
5.3 根据公式 (8) 计算 u(n) 和 v(n) 。
5.4判断是否满足容差条件 (10) ，如果满足，则迭代结束，u(n) 和 v(n)为所求目标值。
5.5 判断是否满足 n <= 最大迭代次数，如果满足，则迭代结束，u(n) 和 v(n)为所求目标值；否则迭代次数+1，转到 3 执行。
6 光流法的灵活性
通过上面演示的求解过程，可知光流法在算法的实现上有很不错的灵活性，主要体现在：
6.1 通过调整 λ 值可以控制计算精度和迭代收敛速度。
6.2 可以只计算图像中的部分子块来计算运动参数。
二、源代码
function varargout = optical_flow(varargin)
% OPTICAL_FLOW M-file for optical_flow.fig
% OPTICAL_FLOW, by itself, creates a new OPTICAL_FLOW or raises the existing
% singleton*.
%
% H = OPTICAL_FLOW returns the handle to a new OPTICAL_FLOW or the handle to
% the existing singleton*.
%
% OPTICAL_FLOW('CALLBACK',hObject,eventData,handles,...) calls the local
% function named CALLBACK in OPTICAL_FLOW.M with the given input arguments.
%
% OPTICAL_FLOW('Property','Value',...) creates a new OPTICAL_FLOW or raises the
% existing singleton*. Starting from the left, property value pairs are
% applied to the GUI before optical_flow_OpeningFunction gets called. An
% unrecognized property name or invalid value makes property application
% stop. All inputs are passed to optical_flow_OpeningFcn via varargin.
%
% *See GUI Options on GUIDE's Tools menu. Choose "GUI allows only one
% instance to run (singleton)".
%
% See also: GUIDE, GUIDATA, GUIHANDLES
% Copyright 2002-2003 The MathWorks, Inc.
% Edit the above text to modify the response to help optical_flow
% Last Modified by GUIDE v2.5 23-Jul-2012 15:10:21
% Begin initialization code - DO NOT EDIT
gui_Singleton = 1;
gui_State = struct('gui_Name', mfilename, ...
'gui_Singleton', gui_Singleton, ...
'gui_OpeningFcn', @optical_flow_OpeningFcn, ...
'gui_OutputFcn', @optical_flow_OutputFcn, ...
'gui_LayoutFcn', [] , ...
'gui_Callback', []);
if nargin && ischar(varargin{1})
gui_State.gui_Callback = str2func(varargin{1});
end
if nargout
[varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:});
else
gui_mainfcn(gui_State, varargin{:});
end
% End initialization code - DO NOT EDIT
% --- Executes just before optical_flow is made visible.
function optical_flow_OpeningFcn(hObject, eventdata, handles, varargin)
% This function has no output args, see OutputFcn.
% hObject handle to figure
% eventdata reserved - to be defined in a future version of MATLAB
% handles structure with handles and user data (see GUIDATA)
% varargin command line arguments to optical_flow (see VARARGIN)
% Choose default command line output for optical_flow
handles.output = hObject;
% Update handles structure
guidata(hObject, handles);
% UIWAIT makes optical_flow wait for user response (see UIRESUME)
% uiwait(handles.figure1);
% --- Outputs from this function are returned to the command line.
function varargout = optical_flow_OutputFcn(hObject, eventdata, handles)
% varargout cell array for returning output args (see VARARGOUT);
% hObject handle to figure
% eventdata reserved - to be defined in a future version of MATLAB
% handles structure with handles and user data (see GUIDATA)
% Get default command line output from handles structure
varargout{1} = handles.output;
% --- Executes on button press in pushbutton_start.
function pushbutton_start_Callback(hObject, eventdata, handles)
% hObject handle to pushbutton_start (see GCBO)
% eventdata reserved - to be defined in a future version of MATLAB
% handles structure with handles and user data (see GUIDATA)
global FilePath1
global FilePath2
global version
if ( isequal(FilePath1,'No Picture') || isequal(FilePath2,'No Picture') )
errordlg('选择图片出错！','MATLAB error');
return
end
switch version
case { 0, 1, 2 }
optic_flow_brox(FilePath1, FilePath2);
case 3
Horn_WJY(FilePath2, FilePath1);
case 4
Lucas_Kanade(FilePath2, FilePath1);
otherwise
errordlg('选择算法出错！','MATLAB error');
end
% --- Executes on button press in pushbutton_pic1.
function pushbutton_pic1_Callback(hObject, eventdata, handles)
% hObject handle to pushbutton_pic1 (see GCBO)
% eventdata reserved - to be defined in a future version of MATLAB
% handles structure with handles and user data (see GUIDATA)
global FilePath1
[FileName,PathName] = uigetfile({'*.jpg';'*.bmp'},'Select the picture');
if ~isequal(FileName,0)
FilePath1 = fullfile(PathName,FileName);
set(handles.edit_lj1,'String',FilePath1);
img=imread(FilePath1);
axes(handles.axes_pic1);
imshow(img);
end
guidata(hObject, handles);
% --- Executes on button press in pushbutton_pic2.
function pushbutton_pic2_Callback(hObject, eventdata, handles)
% hObject handle to pushbutton_pic2 (see GCBO)
% eventdata reserved - to be defined in a future version of MATLAB
% handles structure with handles and user data (see GUIDATA)
global FilePath2
[FileName,PathName] = uigetfile({'*.jpg';'*.bmp'},'Select the picture');
if ~isequal(FileName,0)
FilePath2 = fullfile(PathName,FileName);
set(handles.edit_lj2,'String',FilePath2);
img=imread(FilePath2);
axes(handles.axes_pic2);
imshow(img);
end
guidata(hObject, handles);
function edit_lj1_Callback(hObject, eventdata, handles)
% hObject handle to edit_lj1 (see GCBO)
% eventdata reserved - to be defined in a future version of MATLAB
% handles structure with handles and user data (see GUIDATA)
% Hints: get(hObject,'String') returns contents of edit_lj1 as text
% str2double(get(hObject,'String')) returns contents of edit_lj1 as a double
% --- Executes during object creation, after setting all properties.
function edit_lj1_CreateFcn(hObject, eventdata, handles)
% hObject handle to edit_lj1 (see GCBO)
% eventdata reserved - to be defined in a future version of MATLAB
% handles empty - handles not created until after all CreateFcns called
% Hint: edit controls usually have a white background on Windows.
% See ISPC and COMPUTER.
if ispc
set(hObject,'BackgroundColor','white');
else
set(hObject,'BackgroundColor',get(0,'defaultUicontrolBackgroundColor'));
end
function edit_lj2_Callback(hObject, eventdata, handles)
% hObject handle to edit_lj2 (see GCBO)
% eventdata reserved - to be defined in a future version of MATLAB
% handles structure with handles and user data (see GUIDATA)
% Hints: get(hObject,'String') returns contents of edit_lj2 as text
% str2double(get(hObject,'String')) returns contents of edit_lj2 as a double
% --- Executes during object creation, after setting all properties.
function edit_lj2_CreateFcn(hObject, eventdata, handles)
% hObject handle to edit_lj2 (see GCBO)
% eventdata reserved - to be defined in a future version of MATLAB
% handles empty - handles not created until after all CreateFcns called
% Hint: edit controls usually have a white background on Windows.
% See ISPC and COMPUTER.
if ispc
set(hObject,'BackgroundColor','white');
else
set(hObject,'BackgroundColor',get(0,'defaultUicontrolBackgroundColor'));
end
function edit_arithmetic_Callback(hObject, eventdata, handles)
% hObject handle to edit_arithmetic (see GCBO)
% eventdata reserved - to be defined in a future version of MATLAB
% handles structure with handles and user data (see GUIDATA)
% Hints: get(hObject,'String') returns contents of edit_arithmetic as text
% str2double(get(hObject,'String')) returns contents of edit_arithmetic as a double
% --- Executes during object creation, after setting all properties.
function edit_arithmetic_CreateFcn(hObject, eventdata, handles)
% hObject handle to edit_arithmetic (see GCBO)
% eventdata reserved - to be defined in a future version of MATLAB
% handles empty - handles not created until after all CreateFcns called
% Hint: edit controls usually have a white background on Windows.
% See ISPC and COMPUTER.
if ispc
set(hObject,'BackgroundColor','white');
else
set(hObject,'BackgroundColor',get(0,'defaultUicontrolBackgroundColor'));
end
三、运行结果
四、备注
版本：2014a
