參考自：https://zhuanlan.zhihu.com/p/205778332
http://pages.cs.wisc.edu/~csverma/CS766_09/Stereo/stereo.html

光度立體

光度立體法，即Photometric Stereo, 最早是由當時在MIT的人工智能實驗室的Robert J. Woodham教授在1978年左右提出。

相關論文:

• 《Photometric stereo: A reflectance map technique for determining surface orientation from image intensity》
• 《Photometric Method for Determining Surface Orientation from Multiple Images》

什么是光度立體

光度立體是一種使用多個光源方向估計表面幾何形狀的方法

這種方法的用途是可以重建出物體表面的法向量，以及物體不同表面點的反射率，最關鍵的是它不像傳統的幾何重建（例如立體匹配）方法那樣需要去考慮圖像的匹配問題，因為所需要做的只是采集三張以上，由不同方向的光照射物體的圖像。這個過程中，物體和相機都不動，因此圖像天然就是對齊的，這使得整個過程非常的簡潔。

為什么使用光度立體

光度立體捕獲反射率數據以進行2D匹配

與其他方法的比較：

光度立體法和傳統幾何重建法的比較：

傳統幾何重建法	光度立體法
准確獲取深度	獲取表面梯度
粗糙表面，梯度方向不同時效果好	連續、光滑的表面效果最好
帶有紋理的表面，具有不同反射的表面效果較好	均勻一致的表面效果最好
需要圖像對齊、困難	不需要圖像對齊
不能獲取表面反射比例	可以獲取表面反射比例
	許多情況下，根據表面方向描述對象行形狀比根據參考屏幕上方的范圍或高度描述更可取

基礎方法

假設：

1. 沒有投影/自陰影或高光
2. 灰度/線性成像
3. 入射光由遠處的單一點光源發出，這樣照射到物體表面每一點的光的方向一致，強度一致
4. 相機的投影是正交投影，圖像上的點的坐標可以直接反映三維物體的表面坐標
5. 表面靜止不動
6. 物體表面具有lambertian反射特性，即它對入射光產生漫反射，在每個方向上反射的光強都是一致的

$$ I\propto cos\theta = N · L $$ 其中：

• \(I ：\)輻射強度
• \(N:\)單位法向量
• \(L:\) 單位光源向量
• \(\theta\)：\(N\)和\(L\)之間的夾角

平面方程：

\[Ax+By+Cz+D = 0 \]

平面法線向量：

\[N = [A,B,C]^T \]

由於\(z = -\frac{A}{C}x - \frac{B}{C}y - \frac{D}{C}\)，所以梯度為：

\[\frac{\partial z}{\partial x} = -\frac{A}{C}\qquad\qquad \frac{\partial z}{\partial y} = - \frac{B}{C} \]

重新縮放表面法線

\[n = [\frac{A}{C},\frac{B}{C},1]^T = [-\frac{\partial z}{\partial x},-\frac{\partial z}{\partial y}, 1]^T \]

上式揭示了法向量與梯度之間的關系

考慮使用一個光源對物體進行成像

朗伯仿射模型：

\[E = \rho L cos\theta \]

其中：

• \(E:\)散射
• \(\rho:\)反射率
• \(L:\)光照強度（可忽略）（注意區分這里的L與上面表示光源向量的L表示的不同含義）

取法向量和光源向量的向量積得到\(cos\theta\)

\[pp_s + qq_s + 1 = \sqrt{p^2+q^2+1}\sqrt{p_s^2+q_s^2+1}cos\theta \]

把上式的\(cos\theta\)代入到朗伯定律

\[E = \rho L \frac{pp_s+qq_s+1}{\sqrt{p^2+q^2+1}\sqrt{p_s^2+q_s^2+1}} \]

如果反射率被重新定義且相機響應是線性的，則

\[I' = \rho'\frac{pp_s+qq_s+1}{\sqrt{p^2+q^2+1}\sqrt{p_s^2+q_s^2+1}} \]

已知量：\(p_s,q_s,I\)

未知量：\(p,q,\rho\)

一個方程求三個未知數，是無法求解的，因此至少需要三個方程，即三個光源去照射物體才可以進行求解。

這里面，\(\rho\)和\(N\)代表了物體表面的特性，而L和l則代表了光源的方向和強度，於是可以分別合並表示如下。

\(N\)和\(L\)都是三維空間中的\(3\times 1\)的單位向量，光源的強度可以用常量1表示。像素值\(t\)可以直接通過圖像獲得，右邊的光源方向可以提前標定求得。那么，就可以在物體和相機都不動的情況下，用至少三個不共面的光線照射物體表面第\(i\)點，從而將\(N\)和 \(\rho\)求解出來：

設三個光源方向向量構成矩陣\(L = [L_1,L_2,L_3]^T\)，對應的三個像素值構成向量\(t_i\)，那么有

\[\rho_i = |L^{-1}t_i|\\ N_i = \frac{1}{\rho_i}L^{-1}t_i \]

由於光源照射物體表面時，有可能會產生陰影，這樣三個光源無法同時照亮的區域就會無法求解出結果。所以通常可以采用更多的光源從不同方向分別照亮物體並成像來解決這個問題。

假設有\(n>=3\)個光源，那么則有：

\[L = [L_1,L_2,L_3,...,L_n]^T\\ t_i = [t_{i1},t_{i2},t_{i3},...,t_{in}]^T = \rho_iLN_i \]

那么

\[\rho_i = |(L^TL)^{-1}L^Tt_i|\\ N_i = \frac{1}{\rho}(L^TL)^{-1}Lt_i \]

幾種不同的情況：

1. $p_s=q_s=0 $相當於光線垂直照射

2. \(p_s=q_s=0.5\) 相當於斜着照射

使用圓錐曲線表示：

使用單位圓表示：

因此，光度立體實際上是使用多個光源方向對物體進行照射。

基本過程

光源方向的標定

先第一步，是進行光源方向的標定，一種可行的方法是在場景中放入一個光滑的球並在不同光源下成像，這樣球的表面高亮處就會反映光源的方向

光源方向\(L\)，球表面法向量\(N\)以及反射方向\(R\)之間的關系如下

有

\[L = 2(N·R)·N - R \]

其中\(R\)是反射方向，取為\((0,0,1)\)，是觀察者到物體的向量。\((P_x,P_y)\)是圖像上最亮點的坐標，\((C_x,C_y)\)是圖像上球心的坐標

那么\(N\)的表達式為：

\[N_x = P_x - C_x\\ N_y = P_y - C_y\\ N_z = \sqrt{R^2 - N^2_x - N^2_y} \]

利用\(N\)和\(R\)，就可以求得光源方向\(L\)

計算法向量圖和深度圖

利用上面講過的原理和公式，可以很容易求得物體表面上每一個被足夠多的光源照亮的點的法向量。接下來的關鍵則是利用法向量圖，獲取深度圖。

投影關系如下圖所示：

物體表面法向量\(N\)與物體表面上的向量\(V_1\)及\(V_2\)垂直，那么有：

\[\begin{eqnarray} V_1 &=& (x+1,y,z_{x+1,y})-(x,y,z_{x,y}) \\&=& (1, 0, z_{x+1,y} - z_{x,y})\\ NV_1 &=& 0\\ (n_x,n_y,n_z) · (1, 0, z_{x+1,y} - z_{x,y}) &=& 0\\ n_x + n_z(z_{x+1,y} - z_{x,y}) &=& 0 \end{eqnarray} \]

以及

\[\begin{eqnarray} V_2 &=& (x,y+1,z_{x,y+1})-(x,y,z_{x,y}) \\&=& (0, 1, z_{x,y+1} - z_{x,y})\\ NV_2 &=& 0\\ (n_x,n_y,n_z) · (0, 1, z_{x,y+1} - z_{x,y}) &=& 0\\ n_y + n_z(z_{x,y+1} - z_{x,y}) &=& 0 \end{eqnarray} \]

上面兩個式子，會構成一個線性表達式：

\[Mz = v \]

其中\(v\)由各個像素點的法向量技術而得，\(M\)則是一個尺寸為\((2\times m, m)\)的稀疏矩陣，其中\(m\)為像素個數，\(z\)則是各個像素點的深度坐標構成

通過使方程對稱且為正定，通過共軛梯度法的最小二乘法來求解這種超定線性方程組，最終得到深度圖

利用光度立體法，可以計算得到反照率圖如下，這里對圖像的每個顏色通道都做了獨立運算，而不是把圖像轉換為灰度圖計算

可以獲取到表面法向量，進一步獲取到相對深度圖，甚至還可以做三維渲染，如下圖所示：