光学
光度学基础知识
![[Pasted image 20250428154805.png|center]]
光度学
计量可见光：研究光的强弱–>人眼（主观）
辐射度学
计量电磁辐射：研究辐射强弱–>功率（物理）
辐射能通量 Radiant energy
计算公式为$\Psi=\int \psi(\lambda)d\lambda$，单位为$W$,$\psi(\lambda)$为辐射能谱密度。
视见函数
光度学中，为了表示人眼对不同波长辐射的敏感度差别，定义了一个函数$V(\lambda)$，称为视见函数，又称为光谱光视效率。明亮环境下，人眼对于波长$\lambda=555nm$（黄绿光）的光最灵敏，因此，将波长$\lambda=555mn$的视见函数定为1。$V(\lambda)$为明亮环境下视见函数。$V’(\lambda)$为光线暗的环境下的视见函数。
![[Pasted image 20250428154841.png|center]]
光通量
光通量是人眼感觉到的辐射功率,人类对不同波长的光感受到的辐射通量不同, 光源的亮度根据人眼对各个波长的敏感度进行加权,称为光通量,等于单位时间内某波段的辐射能量和该波段的相对视见率的乘积。标准单位是流明（lumen，简记为lm）。发光强度为1cd（坎德拉）的各向同性的光源在1sr（立体弧度）的立体角内放射的光通量为1lm。光通量就是用来表示辐射功率经过人眼的视见函数影响后的光谱辐射功率大小的物理量。计算公式为：
$$\tag1 \Phi=K_M \int V(\lambda) \psi(\lambda)d\lambda$$
其中，$K_M=683lm/W$为最大光谱当量（$555nm$）。
立体角 Solid angle
立体角，常用字母$\Omega$表示，是一个物体对特定点的三维空间的角度。它描述的是站在某一点的观察者测量到的物体大小的尺度。以观测点为球心，构造一个单位球面，任意物体投影到该单位球面上的投影面积，即为该物体相对于该观测点的立体角。因此，立体角是单位球面上的一块面积，这和“平面角是单位圆上的一段弧长”类似。
锥体的立体角大小定义为，以锥体的顶点为球心作球面，该锥体在球表面截取的面积与球半径平方之比，单位为球面度（$sr$）。
如图所示，物体表面有无限小的面积$dA$，将$dA$投影到图中蓝色区域$dA’$，该区域对应的立体角$d\omega$的计算公式为 $d\omega=\frac{dA’}{r^2}=\frac{dAcos\theta}{r^2}$.
![[Pasted image 20250428154908.png|center]]
微分立体角 Differential Solid Angles
![[file-20250508211618267.png]]
$$dA=(rd\theta)(r\sin\theta d\phi)=r^2\sin\theta d\theta d\phi$$
$$d\omega=\frac{dA}{r^2}=\sin \theta d \theta d \phi$$
辐射强度 Radiant intensity
辐射强度(Radiant Intensity)是点光源单位立体角内发射的辐射通量,记作 $I_e$ , $I_e=\frac{d\Phi_{e}}{d\omega}$,其中$d\Phi_{e}$是立体角$d\omega$内的辐射通量,辐射强度的单位为$W/sr$。
辐射亮度 Radiance
辐射亮度(Radiance)是面辐射源在单位投影面积上、单位立体角内的辐射通量（面元沿某方向单位投影面积和单位立体角的辐射通量），单位是$W\cdot m^{-2} \cdot sr^{-1}$，单位$nit$，尼特。
$$L=\frac{d\Phi}{d\omega dscos\theta}=\frac{dI}{ds cos\theta} \tag2$$
$$L(p,\omega)=\frac{d^2\Phi(p,\omega)}{d\omega dA\cos \theta}$$
辐射照度 Irradiance
power per unit area
辐射照度(Irradiance)是单位面积上的辐射通量,单位是$W/m^2$。辐射照度的计算公式为：
$$E=\frac{d\Phi}{dA} \tag3$$
其中$d\Phi$为辐射通量，$dA$为对应面积。如图所示，有一辐射强度为$I$的点光源，受照面积的辐射照度带入公式得：
$$E=\frac{Id\omega}{dA}=\frac{I\frac{dAcos\theta}{r^2}}{dA}=\frac{Icos\theta}{r^2} \tag4$$
该公式说明 1)在其他条件不变时,当光源处于表面正上方,即$\theta$为0时,辐射照度最大;2)辐射照度与光源和表面的距离的平方成反比,距离越远,表面辐射照度越小。
实际情况中的光源有一定面积,考虑扩展光源如右图所示,扩展光源的辐射亮度为$L$,面积大小为$dB$,发射角度为$\theta_e$,投射到表面$dA$的光能为$(LdBcos\theta_e)d\omega[W]$, 表面对应的辐射照度的计算公式如下:
$$E=\frac{(LdBcos\theta_e)d\omega}{dA}=L(\frac{dBcos\theta_e}{r^2})cos\theta=Ld\omega_icos\theta [W/m^2] \tag5$$
式中$d\omega_i=\frac{dBcos\theta_e}{r^2}$是表面$dA$观察到的光源的立体角。等式表明接收的辐射照度$E$与光源的辐射亮度$L$成正比。式中看上去增大光源与表面的距离,保持立体角相同时,表面接受的辐照度相同,与距离无关,但实际上更远的距离在相同的立体角中会 观察到更大的光源面积,增大的距离和在同一立体角中对应增大的光源面积才保持 辐射照度不变。
![[Pasted image 20250428155026.png|center]]
场景辐射亮度与图像辐射照度的关系
![[Pasted image 20250428155040.png#pic_center|center|650]]
如图所示，从左至右依次是图像平面、镜头和物体表面,设镜头与像平面之间的有效距离为$f$,即成像系统的有效焦距,物体表面的面元$dA_s$通过镜头汇聚后投射到像平面得到对应面元$dA_{i}$,假设物体表面面元$dA_s$处的法线方向为$\pmb n$,相对视线的角度为$\theta$,视线相对光轴的角度为$\alpha$,物体和透镜的距离为$z$,透镜的直径为$d$。由几何关系可知,物体表面面元$dA_s$相对于镜头的立体角与像平面对应面元$dA_i$相对于镜头的立体角相等,$d\omega_s=d\omega_i$,也即$\frac{dA_{i}cos\alpha}{(f/cos\alpha)^2}=\frac{dA_{s}cos\theta}{(z/cos\alpha)^2}$,整理得到
$$\frac{dA_s}{dA_i}=\frac{cos\alpha}{cos\theta}(\frac{z}{f})^2 \tag6$$
透镜相对于表面面元$dA_s$的立体角计算过程如下:
$$d\omega_L = \frac{\frac{\pi d^2}{4}cos\alpha}{(z/cos\alpha)^2}=\frac{\pi d^2 cos^3\alpha}{4z^2} \tag7$$
设表面面元$dA_s$入射到透镜的辐射通量为$d\Phi$,则根据式(5)的前部分有
$$d\Phi = LdA_scos\theta d\omega_L \tag8$$
根据能量守恒,透镜从表面$dA_s$接受的光通量等于透镜投射到$dA_i$的光通量,图像的辐射照度为:
$$E=\frac{d\Phi}{dA_i} \tag9$$
联立（6）、（7）、（8）和（9）得：
$$E=L\frac{\pi}{4}(\frac{d}{f})^{2}cos^{4}\alpha \tag{10}$$
由上式可得图像照度$E$与场景亮度$L$成正比,与有效 f 数$(\frac{f}{d})$的倒数的平方成正比(也即和镜头面积成正比,和有效焦距的平方成反比),图像的照度以$cos^{4}\alpha$的倍率由中心向边缘衰减。
双向反射分布函数
当表面从光源接受光能时,部分光会被反射。出射辐射亮度与入射辐射照度之比 定义为双向反射分布函数(Bidirectional Reflectance Distribution Function,BRDF)。如图所示，光源相对场景的方向为$(\theta_i ,\phi_I)$,辐射照度为$E(\theta_i ,\phi_i)$,在$(\theta_r , \phi_r)$方向上观察到的辐射亮度为$L(\theta_r ,\phi_r)$,则 BRDF 为:
$$f_r(p,\omega_i \to \omega_r)=f(\theta_i,\phi_i,\theta_r,\phi_r)=\frac{L(\theta_r,\phi_r)}{E(\theta_i,\phi_i)} [sr^{-1}] \tag{11}$$
![[Pasted image 20250428155057.png|center]]
双向反射分布函数描述了场景中各个点的反射率属性,具有非负性和互易性,交换入射和出射方向函数值不变,即$f(\theta_i,\phi_i,\theta_r,\phi_r)=f(\theta_r,\phi_r,\theta_i,\phi_i)$。不同表面的 BRDF 不同。光照射到物体表面向各个方向反射的现象为漫反射,理想漫反射体是朗伯体, 朗伯表面从所有方向看上去一样亮,无论从哪个方向观察都会有相同的表面辐射,朗伯表面的 BRDF 是一个常数,如下式所示:
$$f_{lambert}=\frac{\rho_d}{\pi} (0\le \rho_d \le 1) \tag{12}$$
其中，$\rho_d=0$表示完全黑色的表面,$\rho_d=1$表示完全白色的表面。
![[Pasted image 20250428155110.png|center]]
假设图 3-5 中的表面为朗伯表面,光源的辐射强度为$J$,光源方向为$\pmb s$,表面法向为$\pmb n$,光源与法线方向的夹角为$\theta$,观察方向为$\pmb v$,由朗伯表面的 BRDF 有$L=\frac{\rho_d}{\pi}E$, 又根据式（4）可以得到$E=\frac{Jcos\theta_i}{r^2}=\frac{J}{r^2}(\pmb n \cdot \pmb s)$,可推出公式
$$L=\frac{\rho_{d}}{\pi} \frac{J}{r^2}(\pmb n \cdot \pmb s)=\frac{\rho_d}{\pi} \frac{J}{r^2}cos\theta \tag{13}$$
由式可知,在光源及其距离不变的情况下,出射辐射亮度随入射方向偏离表面法线方向而减小。
反射方程
![[file-20250509150026968.png]]
$$L_r(p,\omega_r)=\int^{}_{H^2}f_r(p,\omega_i \to \omega_r)L_i(p,\omega_i)\cos \theta_i \mathrm{d}w_i$$
光度立体法
光度立体是一种从图像强度值恢复表面三维形状的方法。光度立体法使用不同方向光源对物体表面进行照射,通过相机获取的图像的强度值估计表面的法线方向,进一步通过表面梯度重构得到表面高度。
如图所示,光源方向为$\pmb s$,物体表面一点的法线方向为$\pmb n$,相机的观察方向为$\pmb v$,令观察方向为坐标系$z$轴的负方向。图像强度可表示为光源、表面法线方向和表面反射的函数,即图像强度$I=F(光源，表面法线方向，表面反射)$,函数中的光源包括它的方向、辐射强度和距离。在光度立体中将图像强度、光源和表面反射 当作已知量,将表面法线当作未知量,即可通过函数关系求得表面法线方向。
![[Pasted image 20250428155122.png|center]]
假设成像系统是正交投影，则物体上的点$(x,y,z)$投影到像平面$(u,v)$有$u=x,v=y$,物体表面的高度可表示为$z=f(x,y)$，则物体表面的法向量为$[ \frac {\partial f(x,y)}{\partial x},\frac{\partial f(x,y)}{\partial y},-1 ]^T$,定义$p=\frac{\partial f(x,y)}{\partial x}$,$q=\frac{\partial f(x,y)}{\partial y}$，则法线可写为$[p,q,-1]^T$，单位法线向量为$\pmb n =\frac{[p,q,-1]^T}{\sqrt{p^2+q^2+1}}$,将$\pmb n$写为$[n_x,n_y,n_z]^T$，则梯度与单位法线之间的关系为$p=-\frac{n_x}{n_z},q=-\frac{n_y}{n_z}$。对于给定光源和表面反射特性，图像强度可表示为梯度的函数：
$$I=R(p,q) \tag{14}$$
假设入射光来自远点光源，物体表面为朗伯表面，结合式（13），图像强度可表示为：
$$I=c\frac{\rho_d}{\pi}\frac{J}{r^2}cos\theta_i=\rho \pmb n \cdot \pmb s \tag{15}$$
式中$c$为常数,与相机参数有关,在成像系统中是已知的,$p_d$为表面反射率,$\frac{J}{r^2}$是与光源有关的参数,由于光源和相机是给定的,不妨令$\rho=c\frac{\rho_d}{\pi}\frac{J}{r^2}$,$\rho$为与表面反射属性相关的参数。令光源方向为$[p_S,q_S,-1]^T$,可得: $$I=\rho \pmb n \cdot \pmb s=\frac{\rho (1+pp_S+qq_S)}{\sqrt{1+p^2+q^2}\sqrt{1+p_S^2+q_S^2}}=R(p,q) \tag{16}$$
假设有一半径为 60 的球,球心位于物体坐标系原点,有三个光源方向分别为$p_{S1}=0.7,q_{S1}=0.3\mathrm{,}p_{S2}=-0.610,q_{S2}=0.456\mathrm{,}p_{S3}=-0.090,q_{S3}=-0.756$,假设物体各点的$\rho$为 1,则三个光源下相机分别得到的图像如图所示。
![[Pasted image 20250428155136.png|center]]
以物体表面的梯度$p、q$分别作为横纵坐标,将不同方向光源得到的亮度在同一图中画出等高线图如图所示。
![[Pasted image 20250428155148.png|center]]
若球上一点的的亮度在三个光源下分别为$I_1=0.9107,I_2=0.4740,I_3=0.0842$, 则该点的梯度为三个光源下图像亮度等高线的交点$(p,q)=(0.8 ,1)$,如下图所示。如图中所示,一个光源下同一亮度对应等高线上的无数个点,两条等高线的交点可能有两个,因此,至少需要三个光源来确定表面梯度。
![[Pasted image 20250428155203.png|center]]
当图像表面强度、光源方向已知时,可以求得表面法线方向和表面反射率相关参数。假设点 $(x,y)$ 处有三个光源$\pmb s_1,\pmb s_2,\pmb s_3$分别照射,得到的图像强度对应为$I_1,I_2,I_3$,则由朗伯表面光度立体模型有$I_{1}= \rho \pmb n \cdot \pmb s_{1},I_{2}= \rho \pmb n \cdot \pmb s_{2},I_{3}= \rho \pmb n \cdot \pmb s_{3}$$,其中 ,$${\pmb n}= \left[ n_{x},n_{y},n_{z} \right]^{T}$ ,$\pmb{s}i=\left[s{ix}, s_{iy}, s_{iz} \right]^T$,可写成矩阵形式
$$\left[ \begin{matrix} I_{1} \ I_{2} \ I_{3} \end{matrix} \right]= \rho \left[ \begin{matrix} s_{1x}&s_{1y}&s_{1z} \ s_{2x}&s_{2y}&s_{2z} \ s_{3x}&s_{3y}&s_{3z} \end{matrix} \right] \pmb n \tag{17}$$
令$I= \left[ I_{1},I_{2},I_{3} \right]^{T}$，$S$为光源方向矩阵,有
$$\pmb I=\rho \pmb S \pmb n \tag{18}$$
可得$\rho \pmb n =\pmb S ^{-1} \pmb I$，只有当三个光源方向不共面时$\pmb S$可逆,满足此条件时可求得反射率相关参数$\rho = \left| \pmb S^{-1}\pmb I \right|$,表面法线方向为
$$\pmb n =\frac{1}{\rho} \pmb S ^{-1} \pmb I \tag{19}$$
[!note]
只有在表面点能同时被三个光源照射到时才能用式(19)计算表面法线方向。
当光源数大于 3 时,可得到更好的结果,假设有$K$个方向的光源,对应公式变化为
$$\left[ \begin{matrix} I_{1} \ I_{2} \ \vdots \ I_{K} \end{matrix} \right]= \rho \left[ \begin{matrix} s_{1x}&s_{1y}&s_{1z} \ s_{2x}&s_{2y}&s_{2z} \ \vdots& \vdots& \vdots \ s_{Kx}&s_{Ky}&s_{Kz} \end{matrix} \right] \pmb n \tag{20} $$
上式简写为$\pmb I = \rho \pmb S \pmb n=\pmb S \pmb N$,用最小二乘法估计有$\pmb S^T \pmb I= \pmb S^T \pmb S \pmb N$，可求得
$$\mathbf{N} = \left(\mathbf{S}^{\texttt{T}}\mathbf{S}\right)^{ - 1} \mathbf{S}^{\texttt{T}}\mathbf{I} \tag{21}$$
其中$\left(\mathbf{S}^{\texttt{T}}\mathbf{S}\right)^{\texttt{-1}}\mathbf{S}^{ \texttt{T}}$是伪逆矩阵。
对于朗伯表面,多光源同时照射可等价为一个光源照射。由式(15)得,假设有两个不同方向光源同时照射朗伯表面,由它们分别单独照射时获取的图像强度$I_{1}=c \frac{ \rho_{d}}{ \pi} \frac{J_{1}}{r_{1}^{2}}\pmb n \cdot \pmb {s}{1}= \frac{ \rho{d}}{ \pi}k_{1}\pmb n \cdot \pmb {s}{1}$,$I{2}=c \frac{ \rho_{d}}{ \pi} \frac{J_{2}}{r_{2}^{2}}\pmb n \cdot \pmb {s}{2}= \frac{ \rho{d}}{ \pi}k_{2}\pmb n \cdot \pmb {s}{2}$，同时照射得到的图像强度为$I= \frac{ \rho{d}}{ \pi} \pmb n \cdot \left( k_{1}\pmb s_{1}+k_{2}\pmb s_{2} \right)= \frac{ \rho_{d}}{ \pi} \widehat{k}\pmb n \cdot \widehat{\pmb s}$，可等效为亮度参数为$\hat{k}$,方向为$\hat{ \pmb s}$的点光源。进一步 归纳可得多个点光源或扩展光源同时照射朗伯表面可等效为一个光源。
光度立体标定方法
查表法
查找表法是一种数据驱动的标定方法,也是视触觉传感器用于建立 RGB 强度与法线方向的常用方法。若两物体有相同的 BRDF,则它们表面同一法线方向在同一光源的照射下将有相同强度。视触觉传感器使用球形铝粉作为反射膜涂料,表面材料均一,对于不同的下压物体,由于摄像头拍摄到的图像是反射膜,BRDF是一致的,因此,查找表法适用于视触觉传感器。
查找表法的思路如图所示,由于半球理论上有各个方向上的梯度,采用一球形物体进行标定。球半径和摄像头拍摄图像像素所对应的距离已知,根据球的几何关系,通过球下压后摄像头获取圆的半径,可得球下压的深度,结合球的半径得到表面梯度,也即得到图像的 RGB 强度与梯度之间的关系,令表的索引为 RGB 通道的强度,条目为对应的梯度值即得到查找表。实际标定时将球下压在视野的中间和四周,同一强度对应的梯度在各个位置求平均得到。
![[Pasted image 20250428155223.png|center]]
查找表法的优点是不用对物体表面提出假设,对于具有不同材料属性、反射率的反射膜不用假设表面反射情况,可以直接通过数据驱动获得查找表,此种方法对材料有较好的包容性。缺点是可移植性差,由于材料不同和反射膜制作工艺导致的膜的不稳定性,不同材料的反射膜或同一材料不同传感器的反射膜的查找表存在差异,且传感器的体积小,LED光源近,手动焊接可能存在误差,使得不同传感器之间的光源情况不易保持一致,所以更换传感器需重新进行标定。
神经网络法
由于视触觉传感器不符合无穷远点光度立体模型,当物体下压在传感器的不同位置(中心和四周)时梯度和强度的对应关系会产生变化。计算传感器不同位置的梯度角度与强度的对应如图所示,图中第一行为小球下压后的图片,并用圆圈标出接触区域圆,第二行为 x 方向(水平向右为正)的梯度转化为角度值所对应的 RGB 强度值,通过观察发现标定球下压在视触觉传感器的不同位置 x 方向上的梯度与 RGB 强度之间的关系有所不同,主要原因为光源是近点光源,物体和光源之间的距离会影响图像强度。
![[Pasted image 20250428155241.png|center]]
针对梯度和强度的关系受位置影响的情况,查找表法采用各个位置求平均的方式来降低误差。也可以通过神经网络的方法,将位置和强度作为输入,梯度作为输出, 用 BP 神经网络或卷积神经网络(CNN)建立强度、位置与梯度间的映射。神经网络法的优点是可以结合点的位置信息和强度估计点的梯度,且不受反射膜 表面材料的影响。缺点是训练时需要大量图片作为数据来源。
梯度场重构高度图算法
对于接触物体产生的彩色图像,通过图像强度获取表面各个点的梯度后,再用表面重构算法计算得到表面的高度。梯度是高度的微分,对梯度积分理论上可以求解高度,但由于获得的梯度具有噪声,直接通过积分获得表面高度可能造成误差的积累, 积分路径的选择不同也可能得到不同的结果,可以通过对不同路径的积分求平均可以减小噪声的影响。
Frankot-Chellappa 算法
Frankot-Chellappa 算法采用最小二乘法的思想,试图找到一个表面,使得估计的表面高度$z(x,y)$与测量的表面梯度之间的误差最小化,即最小化下式
$$D= \iint \left( \frac{ \partial z}{ \partial x}-p \right)^{2}+ \left( \frac{ \partial z}{ \partial y}-q \right)^{2}dxdy \tag{22}$$
令$Z(u,v),P(u,v)$和$Q(u,v)$分别为$z(x, y),p(x,y)$和$q(x,y)$的傅里叶变换，公式如下所示
$$Z \left( x,y \right)= \iint_{- \infty}^{+ \infty}z \left( u,v \right) e^{j2 \pi \left( ux+vy \right)}dudv \newline P \left( x,y \right)= \iint_{- \infty}^{+ \infty}p \left( u,v \right) e^{j2 \pi \left( ux+vy \right)}dudv\newline Q \left( x,y \right)= \iint_{- \infty}^{+ \infty}q \left( u,v \right) e^{j2 \pi \left( ux+vy \right)}dudv \tag{23}$$
在傅里叶域中最小化目标函数,即找到$Z(u, v)$使$D$最小化,有$\frac{\partial D}{\partial z}=0$,令
$$E \left( z,p,q,z_{x},z_{y} \right)= \left( \frac{ \partial z}{ \partial x}-p \right)^{2}+ \left( \frac{ \partial z}{ \partial y}-q \right)^{2} \tag{24}$$
根据欧拉-拉格朗日方程得
$$\frac{ \partial E}{ \partial z}- \frac{d}{dx} \frac{ \partial E}{ \partial z_{x}}- \frac{d}{dy} \frac{ \partial E}{ \partial z_{y}}=0 \tag{25}$$
$$\frac{\partial E}{\partial z} =0$$
$$\frac{\partial E}{\partial z_x} = 2(z_x-p)$$
$$\frac{\partial E}{\partial z_x}=2(z_x-q)$$
$$\frac{d}{dx} \frac{\partial E}{\partial z_x} = 2 \frac{d}{dx}(z_X-p)=2(z_{xx}-p_x)$$
$$\frac{d}{dy} \frac{\partial E}{\partial z_y} = 2(z_{yy}-q_y)$$
带入式（25）得
$$0-2(z_{xx}-p_x)-2(z_{yy}-q_y)=0$$
$$z_{xx}-z_{yy}=p_x+q_y$$
即 $$\nabla^2z = \frac{\partial p}{\partial x} +\frac{\partial q}{\partial y}$$
$\nabla^2$为拉普拉斯（Laplace）算子
整理得到泊松方程
$$\nabla^2z= \frac{\partial p}{\partial x} + \frac{\partial q}{\partial y} \tag{26}$$
将式（23）带入（26）得
$$\widehat{Z} \left( u,v \right)= \frac{uP \left( u,v \right)+vQ \left( u,v \right)}{2 \pi j \left( u^{2}+v^{2} \right)} \tag{27}$$ 再对$\widehat Z(u,v)$进行傅里叶逆变换得到估计表面高度$\hat z(x,y)$,实现从梯度场到高度场的重构。
光度立体法的三个假设
- 摄像机的投影是正交的，图像上的点的坐标可以直接反映三维物体的表面坐标。
- 每个光源发出的光束平行且均匀，这意味着必须使用强度均匀的远心照明光源，或者应该使用远点光源;并且光源的强度是已知的。
- 最终的物体必须具有朗伯反射特性，即必须以漫反射的方式反射入射光;并且反射涂层的反射率是已知的。
[!note]
要先对相机的相关参数进行标定