多视几何：对极几何的代数表示--基本矩阵F

对极几何的基本概念中是对对极几何的基本形式进行了描述，但并没有从数学角度对其进行描述，而基本矩阵正是对对极几何的代数描述

1.总述

对极几何描述的就是点x和它的对极线l’之间的关系：l′=Fx，其中，矩阵F称为基本矩阵

下面，分别从几何角度和代数角度对上式进行推导

2.几何推导

从几何角度推导关于对极几何方程，如下图所示

• 图像点x是空间点X在图像1中的投影点，那么，有
  
  x=H1X⇒(1)
• 同理，图像点x’是空间点X在图像1中的投影点，那么，有
  
  x′=H2X⇒(2)
• 根据式（1）可以得到
  
  X=H−11x⇒(3)

• 将式（3）代入式（2）得到
  

  x′=H2H−11x→令H2H−11=HπHπx⇒(4)

• 图像点x对应的基线l’是通过点x’和点e’的直线，那么，有
  

  l′=e′×x′⇒(5)
  
  其中的×为叉乘；

• 将式（4）代入式（5）可以得到
  

  l′=e′×Hπx→令e′×Hπ=FFx

结论（来自多视几何教材）：基本矩阵F可以记为F=[e′]×Hπ，其中，Hπ是从一幅图像到另一幅图像通过任意平面π的转移映射，而且，因为[e′]×的秩为2，Hπ的秩为3，所以，F的秩为3

3.代数推导

仍旧以下图为例，假设两幅图像对应的摄像机矩阵分别为P和P’，光心分别为C和C’

图像点x反投影射线的参数形式为X(λ)=P+x+λC

• λ=0时，代表该射线上一个特殊的点P+x
• λ=+∞时，代表该射线上一个特殊的点C

在摄像机2的作用下，这两个特殊的点分别被投影到图像2中的点x’和e’，即

• x′=P′P+x
• e′=P′C

从而，图像1中的点x的对极线即为

l′=e′×x′=P′C×P′P+x=[P′C]×P′P+x

记[P′C]×P′P+=F，则有l′=Fx

4. 对比几何形式和代数形式的基本矩阵

几何形式

F=[e′]×Hπ

代数形式

F=[P′C]×P′P+

5. 基本矩阵可视化

Zisserman的多视几何教程配套代码中，给出了可视化两视点基本矩阵的交互界面

• vgg_ui中的vgg_gui_F.m为可视化基本矩阵的函数
• vgg_examples中view_fund_ex.m为调用vgg_gui_F.m的代码段

原始图像：

下面是几幅运行结果：

view_fund_ex.m内容如下

%%使用vgg_gui_F例子%%% 读入两幅教堂图像im1 = imread('chapel00.png');im2 = imread('chapel01.png');% 读入基本矩阵F的值（好吧，暂时还不知道如何得来）F = load('chapel.00.01.F');% 在GUI中可是化基本矩阵, 利用鼠标点击图像中某点% 可以在另一幅图像中显示该点对应的极线vgg_gui_F(im1, im2, F')

内容如下

function fig=vgg_gui_F(i1,i2,F)%%   fig=vgg_gui_F(i1,i2,F)%%% Visualizes the fundamental matrix of two views%%IN:%   i1 - Matlab image%   i2 - Matlab image%   F - Fundamental matrix (p1'*F*p2=0). Assumes that image coordiantes%       are 1..width where pixel centers are at integer locations. %%OUT:%   fig - handle to the figureif nargin==3   action='start';else   action=i1;   ud = get(gcf, 'UserData');endif strcmp(action,'start'),   if nargin ~= 3      error('Must give 3 arguments... read the docs.\n');   end   h0 = figure('Color',[0.8 0.8 0.8], ...           'NumberTitle','off', ...           'Name','Play With Fundamental Matrix', ...           'ButtonDownFcn', 'disp(''Click on images'')',...           'WindowButtonUpFcn', 'vgg_gui_F(''none'');',...           'WindowButtonMotionFcn', 'vgg_gui_F(''move'')', ...           'Pointer', 'crosshair', ...           'DoubleBuffer', 'on',...           'Units','normalized');   [pointerShape, pointerHotSpot] = CreatePointer;   set(h0, 'Pointer', 'custom', ...       'PointerShapeCData', pointerShape, ...       'PointerShapeHotSpot', pointerHotSpot);   m0=uimenu('Label', '&Color');   uimenu(m0, 'Label', 'blac&K', 'ForegroundColor', [0 0 0], ...      'Accelerator', 'k', 'Callback', 'vgg_gui_F(''ck'');');   uimenu(m0, 'Label', '&Red', 'ForegroundColor', [1 0 0], ...      'Accelerator', 'r', 'Callback', 'vgg_gui_F(''cr'');');   uimenu(m0, 'Label', '&Green', 'ForegroundColor', [0 1 0], ...      'Accelerator', 'g', 'Callback', 'vgg_gui_F(''cg'');');   uimenu(m0, 'Label', '&Blue', 'ForegroundColor', [0 0 1], ...      'Accelerator', 'b', 'Callback', 'vgg_gui_F(''cb'');');   m1=uimenu('Label', '&Size');   uimenu(m1, 'Label', '&Increase', 'Callback', 'vgg_gui_F(''s+'');', 'Accelerator', '+');   uimenu(m1, 'Label', '&Decrease', 'Callback', 'vgg_gui_F(''s-'');', 'Accelerator', '-');   uimenu(m1, 'Label', '&1', 'Callback', 'vgg_gui_F(''s1'');', 'Accelerator', '1');   uimenu(m1, 'Label', '&2', 'Callback', 'vgg_gui_F(''s2'');', 'Accelerator', '2');   uimenu(m1, 'Label', '&3', 'Callback', 'vgg_gui_F(''s3'');', 'Accelerator', '3');   uimenu(m1, 'Label', '&4', 'Callback', 'vgg_gui_F(''s4'');', 'Accelerator', '4');   uimenu(m1, 'Label', '&5', 'Callback', 'vgg_gui_F(''s5'');', 'Accelerator', '5');   uimenu(m1, 'Label', '&6', 'Callback', 'vgg_gui_F(''s6'');', 'Accelerator', '6');   uimenu(m1, 'Label', '&7', 'Callback', 'vgg_gui_F(''s7'');', 'Accelerator', '7');   uimenu(m1, 'Label', '&8', 'Callback', 'vgg_gui_F(''s8'');', 'Accelerator', '8');   uimenu(m1, 'Label', '&9', 'Callback', 'vgg_gui_F(''s9'');', 'Accelerator', '9');   ah1 = axes('Parent', h0, ...          'Position',[0 0 .5 1]);   h1=imshow(i1); hold on; title('Image 1');   set(h1, 'ButtonDownFcn','vgg_gui_F(''b1'');');   ah2 = axes('Parent',h0, ...          'Position',[.5 0 .5 1], ...          'Tag','Axes2');   h2=imshow(i2); hold on; title('Image 2');   set(h2, 'ButtonDownFcn','vgg_gui_F(''b2'');');   point=plot(-1000, -1000,'EraseMode','xor');   l=plot([-1000, -1001], [-1000 -1000], 'r-','EraseMode','xor');   s1=size(i1); s2=size(i2);   t(:,:,1)=F';  t(:,:,2)=F;  F=t;   ud=struct('h0', h0, 'h',[h1 h2], 'ah', [ah1, ah2], ...         'sizes', [s1(1:2); s2(1:2)], ...         'current', -1, 'color', 'k', 'size', 1, ...         'p', point, 'F', F, 'l', l );   set(h0,'UserData',ud);   if nargout > 0, fig = h0; endelseif strcmp(action, 'move')   if ud.current<0 return; end;   pt=get(ud.ah(ud.current),'CurrentPoint');   pts=get_line_points(ud.F(:,:,ud.current)*pt(1,:)', ...               ud.sizes(ud.current,:));   set(ud.l, 'XData', pts(1,:), 'YData', pts(2,:))   set(ud.p, 'XData', pt(1,1), 'YData', pt(1,2))elseif action(1)=='b'   if action(2)=='1' ud.current=1;   elseif action(2)=='2' ud.current=2;   else return;   end   pt=get(ud.ah(ud.current),'CurrentPoint');   pts=get_line_points(ud.F(:,:,ud.current)*pt(1,:)', ud.sizes(1,:));   delete(ud.p);   delete(ud.l);   axes(ud.ah(ud.current));   ud.p=plot(pt(1,1), pt(1,2), [ud.color '+'], ...         'MarkerSize', 8+2*ud.size, 'LineWidth', ud.size,...         'EraseMode','xor');   axes(ud.ah(3-ud.current));   ud.l=plot(pts(1,:), pts(2,:), [ud.color '-'], ...         'LineWidth', ud.size, 'EraseMode','xor');elseif action(1)=='c'   ud.color=action(2);   %get(ud.l)   set(ud.l, 'Color', ud.color);   set(ud.p, 'Color', ud.color);elseif action(1)=='s'   if action(2)=='+' ud.size=ud.size+1;   elseif action(2)=='-' ud.size=max(1, ud.size-1);   else      ud.size = str2num(action(2));   end   set(ud.p, 'LineWidth', ud.size, 'MarkerSize', 8+2*ud.size);   set(ud.l, 'LineWidth', ud.size, 'MarkerSize', 8+2*ud.size);elseif strcmp(action, 'none')   ud.current = -1;else   error(['Unknown command: ' action]);endset(ud.h0, 'UserData',ud);function pts=get_line_points(l,sz)a=l(1); b=l(2);c=l(3);h=sz(1); w=sz(2);% This might cause 'divide by zero' warning:ys=c/-b ;yf=-(a*w+c)/b;xs=c/-a;xf=-(b*h+c)/a;m1 = [[xs;1] [xf;h] [1;ys] [w;yf]];w2 = [(xs<=w & xs>=1) (xf<=w & xf>=1) (ys<=h & ys>=1) (yf<=h & yf>=1)];v = w2>0;pts = [m1(:,v)];% Taken from pixval.m:function [pointerShape, pointerHotSpot] = CreatePointerpointerHotSpot = [8 8];pointerShape = [ ...   NaN NaN NaN NaN NaN   1   2  NaN  2   1 NaN NaN NaN NaN NaN NaN   NaN NaN NaN NaN NaN   1   2  NaN  2   1 NaN NaN NaN NaN NaN NaN   NaN NaN NaN NaN NaN   1   2  NaN  2   1 NaN NaN NaN NaN NaN NaN   NaN NaN NaN NaN NaN   1   2  NaN  2   1 NaN NaN NaN NaN NaN NaN   NaN NaN NaN NaN NaN   1   2  NaN  2   1 NaN NaN NaN NaN NaN NaN    1   1   1   1   1    1   2  NaN  2   1   1   1   1   1   1   1    2   2   2   2   2    2  NaN NaN NaN  2   2   2   2   2   2   2   NaN NaN NaN NaN NaN  NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN    2   2   2   2   2    2  NaN NaN NaN  2   2   2   2   2   2   2    1   1   1   1   1    1   2  NaN  2   1   1   1   1   1   1   1   NaN NaN NaN NaN NaN   1   2  NaN  2   1 NaN NaN NaN NaN NaN NaN   NaN NaN NaN NaN NaN   1   2  NaN  2   1 NaN NaN NaN NaN NaN NaN   NaN NaN NaN NaN NaN   1   2  NaN  2   1 NaN NaN NaN NaN NaN NaN   NaN NaN NaN NaN NaN   1   2  NaN  2   1 NaN NaN NaN NaN NaN NaN   NaN NaN NaN NaN NaN   1   2  NaN  2   1 NaN NaN NaN NaN NaN NaN   NaN NaN NaN NaN NaN NaN NaN  NaN NaN NaN NaN NaN NaN NaN NaN NaN];

---

Reference
[1]《计算机视觉中的多视几何》 ·第八章·对极几何和基本矩阵
[2] 《Multiple View Geometry in Computer Vision Second Edition》·Andrew Zisserman·P239·9 Epipolar Geometry and the Fundamental Matrix
[3] 多视几何教程代码：http://www.robots.ox.ac.uk/~vgg/hzbook/code/