第 4 章　照相机模型与增强现实

第 4 章　照相机模型与增强现实

本章中，我们将尝试对照相机进行建模，并有效地使用这些模型。在之前的章节里，我们已经讲述了图像到图像之间的映射和变换。为了处理三维图像和平面图像之间的映射，我们需要在映射中加入部分照相机产生图像过程的投影特性。本章中我们将会讲述如何确定照相机的参数，以及在具体应用中，如增强现实，如何使用图像间的投影变换。下一章中，我们将使用照相机模型处理其他一些应用，比如多视图及其映射。

4.1　针孔照相机模型

针孔照相机模型（有时称为射影照相机模型）是计算机视觉中广泛使用的照相机模型。对于大多数应用来说，针孔照相机模型简单，并且具有足够的精确度。这个名字来源于一种类似暗箱机的照相机。该照相机从一个小孔采集射到暗箱内部的光线。在针孔照相机模型中，在光线投影到图像平面之前，从唯一一个点经过，也就是照相机中心 C。图 4-1 为从照相机中心前画出图像平面的图解。事实上，在真实的照相机中，图像平面位于照相机中心之后，但是照相机的模型和图 4-1 的模型是一样的。

第 4 章　照相机模型与增强现实 - 图1

图 4-1：针孔照相机模型。图像点 x 是由图像平面与连接三维点 X 和照相机中心 C 的直线相交而成的。虚线表示该照相机的光学坐标轴

由图像坐标轴和三维坐标系中的 x 轴和 y 轴对齐平行的假设，我们可以得出针孔照相机的投影性质。照相机的光学坐标轴和 z 轴一致，该投影几何可以简化成相似三角形。在投影之前通过旋转和平移变换，对该坐标系加入三维点，会出现完整的投影变换。感兴趣的读者可以查阅文献 [13]、[25] 和 [26]。

在针孔照相机中，三维点 X 投影为图像点 x（两个点都是用齐次坐标表示的），如下所示：

λx = P X　　　　　　（4.1）

这里，3×4 的矩阵 P 为照相机矩阵（或投影矩阵）。注意，在齐次坐标系中，三维点 X 的坐标由 4 个元素组成，X=[X, Y, Z, W]。这里的标量 λ 是三维点的逆深度。如果我们打算在齐次坐标中将最后一个数值归一化为 1，那么就会使用到它。

4.1.1　照相机矩阵

照相机矩阵可以分解为：

P = K[R|t]　　　　　　（4.2）

其中，R 是描述照相机方向的旋转矩阵，t 是描述照相机中心位置的三维平移向量，内标定矩阵 K 描述照相机的投影性质。

标定矩阵仅和照相机自身的情况相关，通常情况下可以写成：

$\boldsymbol{K}=\begin{bmatrix}\alpha f & s & c_x\\0 & f & c_y\\0 & 0 & 1\end{bmatrix}$

图像平面和照相机中心间的距离为焦距 f。当像素数组在传感器上偏斜的时候，需要用到倾斜参数 s。在大多数情况下，s 可以设置成 0。也就是说：

$\boldsymbol{K}=\begin{bmatrix}f_x & 0 & c_x\\0 & f_y & c_y\\0 & 0 & 1\end{bmatrix}$ 　　　　　　（4.3）

这里，我们使用了另外的记号 f_x和 f_y，两者关系为 f_x=αf_y。

纵横比例参数 α 是在像素元素非正方形的情况下使用的。通常情况下，我们可以默认设置 α=1。经过这些假设，标定矩阵变为：

$\boldsymbol{K}=\begin{bmatrix}f & 0 & c_x\\0 & f & c_y\\0 & 0 & 1\end{bmatrix}$

除焦距之外，标定矩阵中剩余的唯一参数为光心（有时称主点）的坐标 c=[c_x，c_y]，也就是光线坐标轴和图像平面的交点。因为光心通常在图像的中心，并且图像的坐标是从左上角开始计算的，所以光心的坐标常接近于图像宽度和高度的一半。特别强调一点，在这个例子中，唯一未知的变量是焦距 f。

4.1.2　三维点的投影

下面来创建照相机类，用来处理我们对照相机和投影建模所需要的全部操作：

from scipy import linalg
class Camera(object):
  """ 表示针孔照相机的类"""
  def __init__(self,P):
    """ 初始化 P = K[R|t] 照相机模型"""
    self.P = P
    self.K = None # 标定矩阵
    self.R = None # 旋转
    self.t = None # 平移
    self.c = None # 照相机中心
  def project(self,X):
    """ X（4×n 的数组）的投影点，并且进行坐标归一化 """
  x = dot(self.P,X)
  for i in range(3):
    x[i] /= x[2]
  return x

下面的例子展示如何将三维中的点投影到图像视图中。在这个例子中，我们将使用牛津多视图数据集中的“Model Housing”数据集，可以从 http://www.robots.ox.ac.uk/~vgg/data/data-mview.html 下载。下载这些三维几何图像文件，然后将 house.p3d 文件复制到你的工作目录里：

import camera
# 载入点
points = loadtxt('house.p3d').T
points = vstack((points,ones(points.shape[1])))
# 设置照相机参数
P = hstack((eye(3),array([[0],[0],[-10]])))
cam = camera.Camera(P)
x = cam.project(points)
# 绘制投影
figure()
plot(x[0],x[1],'k.')
show()

首先，我们使用齐次坐标来表示这些点。然后我们使用一个投影矩阵来创建 Camera 对象将这些三维点投影到图像平面并执行绘制操作，输出结果如图 4-2 中间图像所示。

第 4 章　照相机模型与增强现实 - 图5

图 4-2：投影三维点示例：样本图像（左）；图像视图中投影后的点（中）；经过照相机旋转后，投影点的轨迹（右）。数据来自于牛津“Model House”数据集

为了研究照相机的移动会如何改变投影的效果，你可以使用下面的代码。该代码围绕一个随机的三维向量，进行增量旋转的投影。

# 创建变换
r = 0.05*random.rand(3)
rot = camera.rotation_matrix(r)
# 旋转矩阵和投影
figure()
for t in range(20):
  cam.P = dot(cam.P,rot)
  x = cam.project(points)
  plot(x[0],x[1],'k.')
  show()

在上面的代码中，我们使用了 rotation_matrix() 函数，该函数能够创建围绕一个向量进行三维旋转的旋转矩阵（将该函数添加到 camera.py 文件中）：

def rotation_matrix(a):
  """ 创建一个用于围绕向量a 轴旋转的三维旋转矩阵 """
  R = eye(4)
  R[:3,:3] = linalg.expm([[0,-a[2],a[1]],[a[2],0,-a[0]],[-a[1],a[0],0]])
  return R

图 4-2 所示分别为序列中的一幅图像、三维点的投影，以及这些点围绕一个随机向量旋转后三维点投影的轨迹。运行该代码几次，进行不同的随机旋转之后，你会对点在投影中如何旋转有一些感觉。

4.1.3　照相机矩阵的分解

如果给定如方程（4.2）所示的照相机矩阵 P，我们需要恢复内参数 K 以及照相机的位置 t 和姿势 R。矩阵分块操作称为因子分解。这里，我们将使用一种矩阵因子分解的方法，称为 RQ 因子分解。

将下面的方法添加到 Camera 类中：

def factor(self):
  """ 将照相机矩阵分解为K、R、t，其中，P = K[R|t] """
  # 分解前3×3 的部分
  K,R = linalg.rq(self.P[:,:3])
  # 将K 的对角线元素设为正值
  T = diag(sign(diag(K)))
  if linalg.det(T) < 0:
    T[1,1] *= -1
  self.K = dot(K,T)
  self.R = dot(T,R) # T 的逆矩阵为其自身
  self.t = dot(linalg.inv(self.K),self.P[:,3])
  return self.K, self.R, self.t

RQ 因子分解的结果并不是唯一的。在该因子分解中，分解的结果存在符号二义性。由于我们需要限制旋转矩阵 R 为正定的（否则，旋转坐标轴即可），所以如果需要，我们可以在求解到的结果中加入变换 T 来改变符号。

在示例照相机上运行下面的代码，观察照相机矩阵分解的效果：

import camera
K = array([[1000,0,500],[0,1000,300],[0,0,1]])
tmp = camera.rotation_matrix([0,0,1])[:3,:3]
Rt = hstack((tmp,array([[50],[40],[30]])))
cam = camera.Camera(dot(K,Rt))
print K,Rt
print cam.factor()

你会在控制台上得到相同的输出。

4.1.4　计算照相机中心

给定照相机投影矩阵 P，我们可以计算出空间上照相机的所在位置。照相机的中心 C，是一个三维点，满足约束 P C=0。对于投影矩阵为 P=K[R|t] 的照相机，有：

K[R|t]C = K RC + Kt = 0

照相机的中心可以由下述式子来计算：

C = - R^T t

注意，如预期一样，照相机的中心和内标定矩阵 K 无关。

下面的代码可以按照上面公式计算照相机的中心。将其添加到 Camera 类中，该方法会返回照相机的中心：

def center(self):
  """ 计算并返回照相机的中心"""
  if self.c is not None:
    return self.c
else:
    # 通过因子分解计算c
    self.factor()
    self.c = -dot(self.R.T,self.t)
    return self.c

上面的一些方法构成了 Camera 类的基本函数操作。现在，让我们一起探讨如何使用针孔照相机模型。

4.2　照相机标定

标定照相机是指计算出该照相机的内参数。在我们的例子中，是指计算矩阵 K。如果你的应用要求高精度，那么可以扩展该照相机模型，使其包含径向畸变和其他条件。对于大多数应用来说，公式（4.3）中的简单照相机模型已经足够。标定照相机的标准方法是，拍摄多幅平面棋盘模式的图像，然后进行处理计算。例如，OpenCV 中的标定工具使用了这种方法（详见文献 [3]）。

4.2.1　一个简单的标定方法

这里我们将要介绍一个简单的照相机标定方法。大多数参数可以使用基本的假设来设定（正方形垂直的像素，光心位于图像中心），比较难处理的是获得正确的焦距。对于这种标定方法，你需要准备一个平面矩形的标定物体（一个书本即可）、用于测量的卷尺和直尺，以及一个平面。下面是具体操作步骤：

测量你选定矩形标定物体的边长 dX和 dY；
将照相机和标定物体放置在平面上，使得照相机的背面和标定物体平行，同时物体位于照相机图像视图的中心，你可能需要调整照相机或者物体来获得良好的对齐效果；
测量标定物体到照相机的距离 dZ；
拍摄一副图像来检验该设置是否正确，即标定物体的边要和图像的行和列对齐；
使用像素数来测量标定物体图像的宽度和高度 dx 和 dy。

实验设置情况如图 4-3 所示。现在，使用下面的相似三角形（参见图 4-1）关系可以获得焦距：

$f_x=\frac{\text{d}x}{\text{d}X}\text{d}Z, f_y=\frac{\text{d}y}{\text{d}Y}\text{d}Z$

第 4 章　照相机模型与增强现实 - 图7

图 4-3：简单的照相机标定设置：进行标定实验使用的设置情况图像（左）；用于标定的图像（右）。在图像中测量标定物体的宽度和高度，以及设置中标定物体的实际物理尺寸，就可以确定焦距的大小

对于如图 4-3 所示的特定设置，物体宽度和高度的测量值分别为 130 mm 和 185 mm，则，dX=130，dY=185。从照相机到物体的距离为 460 mm，则 dZ=460。你可以使用任意的测量单位，只有测量值的比例才影响最终焦距的计算。你可以使用 ginput() 函数来获得图像中的 4 个点，图像中物体的宽度和高度分别为 722 和 1040 像素。将这些值代入上面的关系表达式可以获得焦距的大小：

f_x= 2555, f_y= 2586

值得注意的是，我们现在获取的焦距是在特定图像分辨率下计算出来的。在这个例子中，图像大小为 2592×1936 像素。记住，焦距和光心是使用像素来度量的，其尺度和图像分辨率相关。如果你使用其他的图像分辨率来拍摄（例如，一个缩略图像），那么这些值都会改变。将照相机的这些测量数值写入到一个如下所示的辅助函数中，会方便一些：

def my_calibration(sz):
  row,col = sz
  fx = 2555*col/2592
  fy = 2586*row/1936
  K = diag([fx,fy,1])
  K[0,2] = 0.5*col
  K[1,2] = 0.5*row
  return K

该函数的输入参数为表示图像大小的元组，返回参数为标定矩阵。这里，我们假设光心为图像的中心。如果你喜欢的话，可以按照焦距的意义来替换焦距的值；对于大多数普通照相机来说，这样做是可以的。注意，这里的标定是对于横向旋转的图像。对于纵向旋转，你需要交换标定矩阵中焦距的值。我们会保留上面的这个函数，方便在下面的章节中使用。

4.3　以平面和标记物进行姿态估计

在第 3 章中，我们学习了如何从平面间估计单应性矩阵。如果图像中包含平面状的标记物体，并且已经对照相机进行了标定，那么我们可以计算出照相机的姿态（旋转和平移）。这里的标记物体可以为对任何平坦的物体。

我们使用一个例子来演示如何进行姿态估计。这里借助图 4-4 中顶端的两幅图像。我们使用下面的代码来提取两幅图像的 SIFT 特征，然后使用 RANSAC 算法稳健地估计单应性矩阵：

import homography
import camera
import sift
# 计算特征
sift.process_image('book_frontal.JPG','im0.sift')
l0,d0 = sift.read_features_from_file('im0.sift')
sift.process_image('book_perspective.JPG','im1.sift')
l1,d1 = sift.read_features_from_file('im1.sift')
# 匹配特征，并计算单应性矩阵
matches = sift.match_twosided(d0,d1)
ndx = matches.nonzero()[0]
fp = homography.make_homog(l0[ndx,:2].T)
ndx2 = [int(matches[i]) for i in ndx]
tp = homography.make_homog(l1[ndx2,:2].T)
model = homography.RansacModel()
H = homography.H_from_ransac(fp,tp,model)

现在我们得到了单应性矩阵。该单应性矩阵将一幅图像中标记物（在这个例子中，标记物是指书本）上的点映射到另一幅图像中的对应点。下面我们定义相应的三维坐标系，使标记物在 X-Y 平面上（Z=0），原点在标记物的某位置上。

为了检验单应性矩阵结果的正确性，我们需要将一些简单的三维物体放置在标记物上，这里我们使用一个立方体。你可以使用下面的函数来产生立方体上的点：

def cube_points(c,wid):
  """ 创建用于绘制立方体的一个点列表（前5 个点是底部的正方形，一些边重合了）"""
  p = []
  # 底部
  p.append([c[0]-wid,c[1]-wid,c[2]-wid])
  p.append([c[0]-wid,c[1]+wid,c[2]-wid])
  p.append([c[0]+wid,c[1]+wid,c[2]-wid])
  p.append([c[0]+wid,c[1]-wid,c[2]-wid])
  p.append([c[0]-wid,c[1]-wid,c[2]-wid]) # 为了绘制闭合图像，和第一个相同
  # 顶部
  p.append([c[0]-wid,c[1]-wid,c[2]+wid])
  p.append([c[0]-wid,c[1]+wid,c[2]+wid])
  p.append([c[0]+wid,c[1]+wid,c[2]+wid])
  p.append([c[0]+wid,c[1]-wid,c[2]+wid])
  p.append([c[0]-wid,c[1]-wid,c[2]+wid]) # 为了绘制闭合图像，和第一个相同
  # 竖直边
  p.append([c[0]-wid,c[1]-wid,c[2]+wid])
  p.append([c[0]-wid,c[1]+wid,c[2]+wid])
  p.append([c[0]-wid,c[1]+wid,c[2]-wid])
  p.append([c[0]+wid,c[1]+wid,c[2]-wid])
  p.append([c[0]+wid,c[1]+wid,c[2]+wid])
  p.append([c[0]+wid,c[1]-wid,c[2]+wid])
  p.append([c[0]+wid,c[1]-wid,c[2]-wid])
  return array(p).T

在上面的函数中，一些数据点会重复出现，plot() 函数会绘制出漂亮的立方体。

第 4 章　照相机模型与增强现实 - 图8

图 4-4：使用平面物体作为标记物，来计算用于新视图投影矩阵的例子。将图像的特征和对齐后的标记匹配，计算出单应性矩阵，然后用于计算照相机的姿态。带有一个灰色正方形区域的模板图像（左上）；从未知视角拍摄的一幅图像，该图像包含同一个正方形，该正方形已经经过估计的单应性矩阵进行了变换（右上）；使用计算出的照相机矩阵变换立方体（下）

有了单应性矩阵和照相机的标定矩阵，我们现在可以得出两个视图间的相对变换：

# 计算照相机标定矩阵
K = my_calibration((747,1000))
# 位于边长为0.2，z=0 平面上的三维点
box = cube_points([0,0,0.1],0.1)
# 投影第一幅图像上底部的正方形
cam1 = camera.Camera( hstack((K,dot(K,array([[0],[0],[-1]])) )) )
# 底部正方形上的点
box_cam1 = cam1.project(homography.make_homog(box[:,:5]))
# 使用H 将点变换到第二幅图像中
box_trans = homography.normalize(dot(H,box_cam1))
# 从cam1 和H 中计算第二个照相机矩阵
cam2 = camera.Camera(dot(H,cam1.P))
A = dot(linalg.inv(K),cam2.P[:,:3])
A = array([A[:,0],A[:,1],cross(A[:,0],A[:,1])]).T
cam2.P[:,:3] = dot(K,A)
# 使用第二个照相机矩阵投影
box_cam2 = cam2.project(homography.make_homog(box))
# 测试：将点投影在z=0 上，应该能够得到相同的点
point = array([1,1,0,1]).T
print homography.normalize(dot(dot(H,cam1.P),point))
print cam2.project(point)

这里我们使用图像的分辨率为 747×1000，第一个产生的标定矩阵就是在该图像分辨率大小下的标定矩阵。下面，我们在原点附近创建立方体上的点。cube_points() 函数产生的前五个点对应于立方体底部的点，在这个例子中对应于位于标记物上 Z=0 平面内的点。第一幅图像（图 4-4 左上）是书本的主视图，我们将其作为这个例子中的模板图像。因为场景坐标的尺度是任意的，所以我们使用下面的矩阵来创建第一个照相机：

$\boldsymbol{P}_1=\boldsymbol{K}\begin{bmatrix}1 & 0 & 0 & 0\\0 & 1 & 0 & 0\\0 & 0 & 1 & -1\end{bmatrix}$

其中，图像的坐标轴和照相机是对齐的，并且放置在标记物的正上方。将前 5 个三维点投影到图像上。有了估计出的单应性矩阵，我们可以将其变换到第二幅图像上。绘

制出变换后的图像，并在同样的标记物位置绘制出这些点（如图 4-4 右上所示）。

现在，结合 P₁ 和 H 构建第二幅图像的照相机矩阵：

P₂ =HP₁

该矩阵将标记平面 Z=0 上的点变换到正确的位置。也就是说，P₂ 矩阵的前两列和第四列是正确的。由于我们知道前 3×3 矩阵块应该为 KR，并且 R 是旋转矩阵，所以我们可以将 P₂ 乘以标定矩阵的逆，然后将第三列替换为前两列的交叉乘积，以此来恢复第三列。

作为合理性验证，我们可以使用新矩阵投影标记平面的一个点，然后检查投影后的点是否与使用第一个照相机和单应性矩阵变换后的点相同。你会发现，控制台上得到了相同的输出结果。

你可以用如下所示的代码来可视化这些投影后的点：

im0 = array(Image.open('book_frontal.JPG'))
im1 = array(Image.open('book_perspective.JPG'))
# 底部正方形的二维投影
figure()
imshow(im0)
plot(box_cam1[0,:],box_cam1[1,:],linewidth=3)
# 使用H 对二维投影进行变换
figure()
imshow(im1)
plot(box_trans[0,:],box_trans[1,:],linewidth=3)
# 三维立方体
figure()
imshow(im1)
plot(box_cam2[0,:],box_cam2[1,:],linewidth=3)
show()

该代码绘制出如图 4-4 所示的三幅图像。为了能够在后面的例子中再次使用这些计算的结果，我们可以使用 Pickle 将这些照相机矩阵保存起来：

import pickle
with open('ar_camera.pkl','w') as f:
  pickle.dump(K,f)
  pickle.dump(dot(linalg.inv(K),cam2.P),f)

现在我们已经学会计算给定平面场景物体的照相机矩阵。我们结合特征匹配和单应性矩阵，以及照相机标定技术，简单演示了如何在一幅图像上放置一个立方体。有了照相机的姿态估计技术，我们现在就具备了创建简单增强现实应用的基本技能了。

4.4　增强现实

增强现实（Augmented Reality，AR）是将物体和相应信息放置在图像数据上的一系列操作的总称。最经典的例子是放置一个三维计算机图形学模型，使其看起来属于该场景；如果在视频中，该模型会随着照相机的运动很自然地移动。如上一节所示，给定一幅带有标记平面的图像，我们能够计算出照相机的位置和姿态，使用这些信息来放置计算机图形学模型，能够正确表示它们。在本章的最后一节，我们将介绍如何建立一个简单的增强现实例子。其中，我们会用到两个工具包：PyGame 和 PyOpenGL。

4.4.1　PyGame和PyOpenGL

PyGame 是非常流行的游戏开发工具包，它可以非常简单地处理显示窗口、输入设备、事件，以及其他内容。PyGame 是开源的，可以从 http://www.pygame.org/ 下载。事实上，它是一个 Python 绑定的SDL游戏引擎。你可以查看附录 A 来获取关于安装的帮助。你还可以查看文献 [21] 来获取关于 PyGame 工具包的更多编程细节。

PyOpenGL 是 OpenGL 图形编程的 Python 绑定接口。OpenGL 可以安装在几乎所有的系统上，并且具有很好的图形性能。OpenGL 具有跨平台性，能够在不同的操作系统之间工作。关于 OpenGL 的更多信息，参见 http://www.opengl.org/。在开始页面（http://www.opengl.org/wiki/Getting_started），有针对初学者的很多资源。PyOpenGL 是开源的，并且易于安装。你可以从附录 A 了解关于 PyOpenGL 更多内容。你同样可以在项目网页 http://pyopengl.sourceforge.net/ 获取更多细节内容。

我们没有涵盖关于 OpenGL 编程的绝大部分内容。相反，我们在这里只讲述重要的内容，例如如何在 OpenGL 中使用照相机矩阵，如何设置一个基本的三维模型。你可以从 PyOpenGL-Demo 工具包（http://pypi.python.org/pypi/PyOpenGL-Demo）获取一些很好的例子和演示。如果你不熟悉 PyOpenGL，那么该工具包是个很好的开始。

我们使用 OpenGL 将一个三维模型放置在一个场景中。为了使用 PyGame 和 PyOpenGL 工具包来完成该应用，需要在脚本的开始部分载入下面的命令：

from OpenGL.GL import 
from OpenGL.GLU import 
import pygame, pygame.image
from pygame.locals import *

你可以看到，这里主要使用 OpenGL 中的两个部分：GL 部分包含所有以“gl”开头的函数，其中包含我们需要的大部分函数；GLU 部分是 OpenGL 的实用函数库，里面包含了一些高层的函数。我们主要使用它来设置照相机投影。pygame 部分用来设置窗口和事件控制；其中 pygame.image 用来载入图像和创建 OpenGL 的纹理，pygame.locals 用来设置 OpenGL 的显示区域。

需要对一个 OpenGL 场景进行两个部分的设置：投影和视图矩阵的建模。让我们一起来学习如何由针孔照相机来创建这些矩阵。

4.4.2　从照相机矩阵到OpenGL格式

OpenGL 使用 4×4 的矩阵来表示变换（包括三维变换和投影）。这和我们使用的 3×4 照相机矩阵略有差别。但是，照相机与场景的变换分成了两个矩阵，GL_PROJECTION 矩阵和 GL_MODELVIEW 矩阵。GL_PROJECTION 矩阵处理图像成像的性质，等价于我们的内标定矩阵 K。GL_MODELVIEW 矩阵处理物体和照相机之间的三维变换关系，对应于我们照相机矩阵中的 R 和 t 部分。一个不同之处是，假设照相机为坐标系的中心，GL_MODELVIEW 矩阵实际上包含了将物体放置在照相机前面的变换。OpenGL 有很多特性，我们会在下面例子中一一讲解。

假设我们已经获得了标定好的照相机，即已知标定矩阵 K，下面的函数可以将照相机参数转换为 OpenGL 中的投影矩阵：

def set_projection_from_camera(K):
  """ 从照相机标定矩阵中获得视图"""
  glMatrixMode(GL_PROJECTION)
  glLoadIdentity()
  fx = K[0,0]
  fy = K[1,1]
  fovy = 2*arctan(0.5*height/fy)*180/pi
  aspect = (width*fy)/(height*fx)
  # 定义近的和远的剪裁平面
  near = 0.1
  far = 100.0
  # 设定透视
  gluPerspective(fovy,aspect,near,far)
  glViewport(0,0,width,height)

我们假设标定矩阵如（4.3）那样简单，光心为图像的中心。第一个函数 glMatrixMode() 将工作矩阵设置为 GL_PROJECTION，接下来的命令会修改这个矩阵 1。然后，glLoadIdentity() 函数将该矩阵设置为单位矩阵，这是重置矩阵的一般操作。然后，我们根据图像的高度、照相机的焦距以及纵横比，计算出视图中的垂直场。OpenGL 的投影同样具有近距离和远距离的裁剪平面来限制场景拍摄的深度范围。我们设置近深度为一个小的数值，使得照相机能够包含最近的物体，而远深度设置为一个大的数值。我们使用 GLU 的实用函数 gluPerspective() 来设置投影矩阵，将整个图像定义为视图部分（也就是显示的部分）。和下面的模拟视图函数相似，你可以使用 glLoadMatrixf() 函数的一个选项来定义一个完全的投影矩阵。当简单版本的标定矩阵不够好时，可以使用完全投影矩阵。

1这是个奇怪的处理方式，但是只需要处理 GL_PROJECTION 和 GL_MODELVIEW 两个矩阵，所以是可行的。

模拟视图矩阵能够表示相对的旋转和平移，该变换将该物体放置在照相机前（效果是照相机在原点上）。模拟视图矩阵是个典型的 4×4 矩阵，如下所示：

$\begin{bmatrix}\boldsymbol{R} & \boldsymbol{t}\\\boldsymbol{0} & \boldsymbol{1}\end{bmatrix}$

其中，R 是旋转矩阵，列向量表示 3 个坐标轴的方向，t 是平移向量。当创建模拟视图矩阵时，旋转矩阵需要包括所有的旋转（物体和坐标系的旋转），可以将单个旋转分量相乘来获得旋转矩阵。

下面的函数实现如何获得移除标定矩阵后的 3×4 针孔照相机矩阵（将 P 和 K^-1 相乘），并创建一个模拟视图：

def set_modelview_from_camera(Rt):
  """ 从照相机姿态中获得模拟视图矩阵"""
  glMatrixMode(GL_MODELVIEW)
  glLoadIdentity()
  # 围绕x 轴将茶壶旋转90 度，使z 轴向上
  Rx = array([[1,0,0],[0,0,-1],[0,1,0]])
  # 获得旋转的最佳逼近
  R = Rt[:,:3]
  U,S,V = linalg.svd(R)
  R = dot(U,V)
  R[0,:] = -R[0,:] # 改变x 轴的符号
  # 获得平移量
  t = Rt[:,3]
  # 获得4×4 的模拟视图矩阵
  M = eye(4)
  M[:3,:3] = dot(R,Rx)
  M[:3,3] = t
  # 转置并压平以获取列序数值
  M = M.T
  m = M.flatten()
  # 将模拟视图矩阵替换为新的矩阵
  glLoadMatrixf(m)

在上面函数中，我们首先切换到 GL_MODELVIEW 矩阵，然后重置该矩阵。然后，由于需要旋转该物体（你将在下面看到），所以我们创建一个 90 度的旋转矩阵。接下来，由于估计照相机矩阵时，可能会有错误或者噪声干扰，所以我们确保照相机矩阵的旋转部分确实是个旋转矩阵。该操作使用 SVD 分解方法，旋转矩阵的最佳逼近可以通过 R=UV^T 来获得。由于 OpenGL 中的坐标系和上面用到的有点不同，所以我们将 x 轴翻转。然后，我们将模拟视图矩阵 M 设置为旋转矩阵的乘积。glLoadMatrixf() 函数通过输入参数为按列排列的 16 个数值数组，来设置模拟视图。将 M 矩阵转置，然后压平并输入 glLoadMatrixf() 函数。

4.4.3　在图像中放置虚拟物体

我们需要做的第一件事是将图像（打算放置虚拟物体的图像）作为背景添加进来。在 OpenGL 中，该操作可以通过创建一个四边形的方式来完成，该四边形为整个视图。完成该操作最简单的方式是绘制出四边形，同时将投影和模拟试图矩阵重置，使得每一维的坐标范围在 -1 到 1 之间。

下面的函数可以载入一幅图像，然后将其转换成一个 OpenGL 纹理，并将该纹理放置在四边形上：

def draw_background(imname):
  """ 使用四边形绘制背景图像"""
  # 载入背景图像（应该是.bmp 格式），转换为OpenGL 纹理
  bg_image = pygame.image.load(imname).convert()
  bg_data = pygame.image.tostring(bg_image,"RGBX",1)
  glMatrixMode(GL_MODELVIEW)
  glLoadIdentity()
  glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT)
  # 绑定纹理
  glEnable(GL_TEXTURE_2D)
  glBindTexture(GL_TEXTURE_2D,glGenTextures(1))
  glTexImage2D(GL_TEXTURE_2D,0,GL_RGBA,width,height,0,GL_RGBA,GL_UNSIGNED_BYTE,bg_data)
  glTexParameterf(GL_TEXTURE_2D,GL_TEXTURE_MAG_fiLTER,GL_NEAREST)
  glTexParameterf(GL_TEXTURE_2D,GL_TEXTURE_MIN_fiLTER,GL_NEAREST)
  # 创建四方形填充整个窗口
  glBegin(GL_QUADS)
  glTexCoord2f(0.0,0.0); glVertex3f(-1.0,-1.0,-1.0)
  glTexCoord2f(1.0,0.0); glVertex3f( 1.0,-1.0,-1.0)
  glTexCoord2f(1.0,1.0); glVertex3f( 1.0, 1.0,-1.0)
  glTexCoord2f(0.0,1.0); glVertex3f(-1.0, 1.0,-1.0)
  glEnd()
  # 清除纹理
  glDeleteTextures(1)

该函数首先使用 PyGame 中的一些函数来载入一幅图像，将其序列化为能够在 PyOpenGL 中使用的原始字符串表示。然后，重置模拟视图，清除颜色和深度缓存。接下来，绑定这个纹理，使其能够在四边形和指定插值中使用它。四边形是在每一维分别为 -1 和 1 的点上定义的。注意，纹理图像的坐标是从 0 到 1。最后，清除该纹理，避免其干扰之后准备绘制的图像。

现在已经准备好将物体放置入场景中。我们将使用“hello world”的计算机图形学例子，Utah 茶壶（http://en.wikipedia.org/wiki/Utah_teapot）。这个茶壶有丰富的历史，在 GLUT 用作一个标准形状：

from OpenGL.GLUT import *
glutSolidTeapot(size)

该命令产生一个相对大小为 size 的茶壶模型。

下面的函数将会设置颜色和其他特性，产生一个红色的漂亮茶壶：

def draw_teapot(size):
  """ 在原点处绘制红色茶壶"""
  glEnable(GL_LIGHTING)
  glEnable(GL_LIGHT0)
  glEnable(GL_DEPTH_TEST)
  glClear(GL_DEPTH_BUFFER_BIT)
  # 绘制红色茶壶
  glMaterialfv(GL_FRONT,GL_AMBIENT,[0,0,0,0])
  glMaterialfv(GL_FRONT,GL_DIFFUSE,[0.5,0.0,0.0,0.0])
  glMaterialfv(GL_FRONT,GL_SPECULAR,[0.7,0.6,0.6,0.0])
  glMaterialf(GL_FRONT,GL_SHININESS,0.25*128.0)
  glutSolidTeapot(size)

上面的函数中，前两行激活了灯光效果和一盏灯。灯被计为 GL_LIGHT0 和 GL_LIGHT1 等。在本例中，我们只使用一盏灯。glEnable() 函数用来激活 OpenGL 的一些特性。这些特性是使用大写常量来定义的。关闭特性可以使用相应的 glDisable() 函数来完成。接下来，深度测试被激活，使物体按照其深度表示出来（远处的物体不能绘制在近处物体的前面），然后清理深度缓存。接下来，指定物体的物质特性，例如漫反射和镜面反射颜色。在最后一行代码中，将指定的物质特性加入到 Utah 茶壶上。

4.4.4　综合集成

下面的完整脚本可以生成如图 4-5 所示的图像（假设你已经将上面的函数添加到了同一个文件中）：

from OpenGL.GL import 
from OpenGL.GLU import 
from OpenGL.GLUT import 
import pygame, pygame.image
from pygame.locals import 
import pickle
width,height = 1000,747
def setup():
  """ 设置窗口和pygame 环境"""
  pygame.init()
  pygame.display.set_mode((width,height),OPENGL | DOUBLEBUF)
  pygame.display.set_caption('OpenGL AR demo')
# 载入照相机数据
with open('ar_camera.pkl','r') as f:
  K = pickle.load(f)
  Rt = pickle.load(f)
setup()
draw_background('book_perspective.bmp')
set_projection_from_camera(K)
set_modelview_from_camera(Rt)
draw_teapot(0.02)
while True:
  event = pygame.event.poll()
  if event.type in (QUIT,KEYDOWN):
  break
pygame.display.flip()

第 4 章　照相机模型与增强现实 - 图11

图 4-5：增强现实。使用由特征匹配计算出的照相机参数，将一个计算机图形学模型放置在场景中的书本上：将 Utah 茶壶按照和坐标轴对齐的方式显示出来（上）；进行合理性验证，查看原点的位置（下）

首先，该脚本使用 Pickle 载入照相机标定矩阵，以及照相机矩阵中的旋转和平移部分。假设这些信息已经按照 4.3 节的描述进行保存。setup() 函数会初始化 PyGame，将窗口设置为图像的大小，绘制图像区域为两倍的 OpenGL 窗口缓存大小。接下来，载入背景图像，使其与窗口相符。然后，设定照相机和模拟视图矩阵。最后，在正确的位置上绘制出茶壶。

在 PyGame 中，使用带有对所有变化进行定期询问的无限循环来处理事件。这些事件可以为键盘、鼠标，或者其他。在这个例子中，我们检查应用是否退出，或者是否有按键按下，如果是，则退出这个循环。pygame.display.flip() 命令将会在屏幕上绘制出物体。

程序的输出结果如图 4-5 所示。可以看到，物体的朝向是正确的（在图 4-4 中，茶壶和立方体的边是对齐的）。为了验证放置的正确性，你可以通过给 size 变量传递一个较小的数值，将茶壶变小。在图 4-4 中，茶壶应该放置在靠近立方体坐标为 [0, 0, 0] 的角上。图 4-5 给出了一个例子。

4.4.5　载入模型

在结束本章之前，我们讲述最后一个细节：载入并显示三维模型。你可以在 http://www.pygame.org/wiki/OBJFileLoader 了解如何在 PyGame 中载入保存在 .obj 格式文件中的模型。除此之外，你还可以在 http://en.wikipedia.org/wiki/Wavefront_.obj_file 学习关于 .obj 以及其他相关文件格式的更多内容。

下面使用一个基本的例子来说明其使用方法。我们将使用一个免费的玩具飞机模型，模型来自 http://www.oyonale.com/modeles.php2。下载其 .obj 文件，然后保存为 toyplane.obj。当然，你也可以用任何其他模型替换，下面的代码不变。

2模型由 Gilles Tran 提供（共享创意许可）。

假设已经下载模型文件并命名为 objloader.py, 将下面的函数添加到茶壶例子的代码文件中：

def load_and_draw_model(filename):
  """ 使用objloader.py，从.obj 文件中装载模型
    假设路径文件夹中存在同名的.mtl 材料设置文件"""
  glEnable(GL_LIGHTING)
  glEnable(GL_LIGHT0)
  glEnable(GL_DEPTH_TEST)
  glClear(GL_DEPTH_BUFFER_BIT)
  # 设置模型颜色
  glMaterialfv(GL_FRONT,GL_AMBIENT,[0,0,0,0])
  glMaterialfv(GL_FRONT,GL_DIFFUSE,[0.5,0.75,1.0,0.0])
  glMaterialf(GL_FRONT,GL_SHININESS,0.25*128.0)
  # 从文件中载入
  import objloader
  obj = objloader.OBJ(filename,swapyz=True)
  glCallList(obj.gl_list)

和前面一样，我们首先设置模型的照明条件和颜色属性。接下来，我们将模型载入一个 obj 对象中，然后在文件中执行 OpenGL 的调用。

你可以在相应的 .mtl 文件中设置纹理和材料属性。实际上，objloader 模块需要一个含有材料设置的文件。我们采用仅创建一个小材料文件的实用方法，而不修改脚本。在这个例子中，我们仅指定了颜色信息。

使用下面的命令来创建 toyplane.mtl 文件：

newmtl lightblue
Kd 0.5 0.75 1.0
illum 1

上面的命令将物体的漫反射颜色设置为浅灰蓝色。现在，确保将 .obj 文件中的“usemtl”标签替换为：

usemtl lightblue

我们将添加纹理留作练习。将上面例子中的 draw_teapot() 命令替换为：

load_and_draw_model('toyplane.obj')

程序就会生成如图 4-6 所示的窗口图像。

这是本书关于增强现实和 OpenGL 最深入的例子。学习了标定矩阵、计算照相机姿态、将照相机转变成 OpenGL 格式，以及在场景中显示模型等方法，你已经具备了继续探索增强现实的基础。在下一节中，我们将在不使用标记物的情况下，继续学习照相机模型，计算三维结构以及照相机姿态。

第 4 章　照相机模型与增强现实 - 图12

图 4-6：从 .obj 文件载入三维模型，并使用由特征匹配计算出的照相机参数将其放置在书本上

练习

修改用于图 4-2 运动的示例代码，使其能够对点而非照相机进行变换，你应该会得到相同的图像。使用不同的变换进行实验，然后绘制出相应的结果。
一些牛津多视图数据集已经给定了照相机矩阵。对一个这样的数据集，计算其照相机位置，并且绘制出照相机的路径。这和你在图像中看到的是否吻合？
拍摄一些带有平面标记物和物体的场景图像。将图像的特征和一幅全景图像的特征相匹配，计算出每幅图像照相机位置的姿态。绘制出照相机的轨迹，以及标记物的平面。如果你喜欢，也可以在图像中加入特征点。
在增强现实的例子中，假设物体放置在原点，并且只将照相机的位置应用到模拟视图矩阵中。将物体之间的变换加入到矩阵中来修改该例子，使得在不同的位置放置一些物体。例如，在标记位置处放置茶壶网格。
浏览 .obj 模型文件的在线文档，学习如何使用纹理模型。找一个模型（或者创建自己的），然后将其添加到场景中。

第 4 章 照相机模型与增强现实

第 4 章 照相机模型与增强现实

4.1 针孔照相机模型

4.1.1 照相机矩阵

4.1.2 三维点的投影

4.1.3 照相机矩阵的分解

4.1.4 计算照相机中心

4.2 照相机标定

4.2.1 一个简单的标定方法

4.3 以平面和标记物进行姿态估计

4.4 增强现实

4.4.1 PyGame和PyOpenGL

4.4.2 从照相机矩阵到OpenGL格式

4.4.3 在图像中放置虚拟物体

4.4.4 综合集成

4.4.5 载入模型

练习

第 4 章　照相机模型与增强现实

第 4 章　照相机模型与增强现实

4.1　针孔照相机模型

4.1.1　照相机矩阵

4.1.2　三维点的投影

4.1.3　照相机矩阵的分解

4.1.4　计算照相机中心

4.2　照相机标定

4.2.1　一个简单的标定方法

4.3　以平面和标记物进行姿态估计

4.4　增强现实

4.4.1　PyGame和PyOpenGL

4.4.2　从照相机矩阵到OpenGL格式

4.4.3　在图像中放置虚拟物体

4.4.4　综合集成

4.4.5　载入模型