标题:基于数字光栅投影的结构光三维测量技术与系统研究 [博]
单位:华中科技大学
作者:李中伟
年份:2009
注:这里只记了相位出现问题的原因分析,如何校正,看最新论文。
相位测量轮廓术,误差主要有四个来源:
-
相移机构的相移误差:
相移步距不相等,在机械式相移中,相移难避免,但DLP数字式光栅这块误差很小。
-
随机噪声:
相移法对随机噪声有较好的抑制,且相位计算前可以通过高斯滤波等图像预处理操作进行滤除。
-
光栅图像条纹,离散化误差:
选择合适的光栅图像周期、灰度范围、投影仪和相机的像素比。
-
光栅图像条纹,非正弦化(主),由两部分引起:投影仪的Gamma非线性曲线、CCD相机的非线性响应
举例:如果投出去的条纹非正弦,那么我们的解码就会有问题,比如说:
- 相邻两个像素点,左边像素投射出去亮度为 148,右边像素同样为 148,那么就很难区分这两个像素点;
- 投射出去条纹非正弦化,解码公式同样存在问题。
01 传递模型
在整个测量过程中,CCD相机拍摄光栅图像的过程如下图所示:
假设由计算机生成并发送到投影仪离线的光栅图像灰度分布:
I
n
(
x
,
y
)
=
a
0
+
b
0
cos
[
φ
(
x
,
y
)
+
δ
n
]
{I_n}\left( {x,y} \right) = {a_0} + {b_0}\cos \left[ {\varphi \left( {x,y} \right) + {\delta _n}} \right]
In(x,y)=a0+b0cos[φ(x,y)+δn]
其中:
- a 0 a_0 a0:平均亮度
- b 0 b_0 b0:调制幅值
- φ ( x , y ) \varphi(x,y) φ(x,y):相位主值
- δ n \delta_n δn:第 i i i 副图像的相移值
经投影仪投射出的光栅灰度分布:
I
n
p
(
x
,
y
)
=
f
p
(
I
n
(
x
,
y
)
)
I_n^p\left( {x,y} \right) = {f_p}\left( {{I_n}(x,y)} \right)
Inp(x,y)=fp(In(x,y))
其中**:
f
p
(
I
)
f_p(I)
fp(I) 为
I
I
I 的函数,代表投影仪对输入光强的响应。**
假设投影仪投射在被测物体上的光栅图像的反射率为
r
(
x
,
y
)
r(x,y)
r(x,y),周围的环境光为
a
1
(
x
,
y
)
a_1(x,y)
a1(x,y),则被测物体反射回来的光的灰度分布:
I
n
o
(
x
,
y
)
=
r
(
x
,
y
)
[
I
n
p
(
x
,
y
)
+
a
1
(
x
,
y
)
]
I_n^o\left( {x,y} \right) = r\left( {x,y} \right)\left[ {I_n^p\left( {x,y} \right) + {a_1}\left( {x,y} \right)} \right]
Ino(x,y)=r(x,y)[Inp(x,y)+a1(x,y)]
注:图案已经变形,灰度分布依旧,所以这里强调的是灰度分布。
经被测物体反射回来的光
I
n
p
(
x
,
y
)
I_n^p(x,y)
Inp(x,y) 被
C
C
D
CCD
CCD 相机拍摄,假设进入相机的环境光为
a
2
(
x
,
y
)
a_2(x,y)
a2(x,y),则相机最终拍摄到的光栅图像的光强分布:
I
n
c
(
x
,
y
)
=
f
c
(
I
n
O
(
x
,
y
)
+
a
2
(
x
,
y
)
)
I_n^c(x,y)=f_c(I_n^O(x,y)+a_2(x,y))
Inc(x,y)=fc(InO(x,y)+a2(x,y))
其中**:
f
c
(
I
)
fc_(I)
fc(I) 为
I
I
I 的函数,代表CCD相机对输入光强的响应。**
联立前面的式子,可得:
I
n
c
(
x
,
y
)
=
f
(
I
n
(
x
,
y
)
)
=
f
c
(
f
p
(
I
n
(
x
,
y
)
)
)
I_n^c(x,y)=f(I_n(x,y))=f_c(f_p(I_n(x,y)))
Inc(x,y)=f(In(x,y))=fc(fp(In(x,y)))
02 误差原因
这里的误差主要有两部分:
- 对于环境光的干扰,我们采用蓝波波段投影,在相机端加相应的滤色片,去除其他波段的光,从而大幅降低环境光的影响
- 投影仪、相机的非线性响应,引起拍摄到的光栅图像的灰度的非正弦化
注:对于物体表面反射率不一致,相移法通过外差的方法,这种影响,在计算相位解码过程中得到了很大程度的抑制。
对于投影仪的非线性响应,其主要由伽马变换引起。在普通商用投影仪中,由于人眼的非线性感光,伽马变换用来改善投影仪投射图案的视觉效果,可以由指数变换表示:
I
n
p
(
x
,
y
)
=
[
I
n
(
x
,
y
)
]
γ
I_n^p(x,y)=[I_n(x,y)]^{\gamma}
Inp(x,y)=[In(x,y)]γ
CCD相机输入光强与输出电压之间是非线性的,影响较大的是二阶、三阶的非线性响应,可以表示:
I
n
c
(
x
,
y
)
=
e
0
+
e
1
[
I
n
0
(
x
,
y
)
+
a
2
(
x
,
y
)
]
+
e
2
[
I
n
0
(
x
,
y
)
+
a
2
(
x
,
y
)
]
2
+
e
3
[
I
n
0
(
x
,
y
)
+
a
2
(
x
,
y
)
]
3
I_n^c\left( {x,y} \right) = {e_0} + {e_1}\left[ {I_n^0\left( {x,y} \right) + {a_2}\left( {x,y} \right)} \right] + {e_2}{\left[ {I_n^0\left( {x,y} \right) + {a_2}\left( {x,y} \right)} \right]^2} + {e_3}{\left[ {I_n^0\left( {x,y} \right) + {a_2}\left( {x,y} \right)} \right]^3}
Inc(x,y)=e0+e1[In0(x,y)+a2(x,y)]+e2[In0(x,y)+a2(x,y)]2+e3[In0(x,y)+a2(x,y)]3
Q:相机为什么出厂前不通过软件做成线性的呢?
可能原因:跟光圈、焦距这些有关,而这些参数是可以调整的,因此无法在出厂就做好。
将前面三个式子结合,滤除掉环境光的干扰,可得:
I
n
c
(
x
,
y
)
=
k
0
+
k
1
[
I
n
(
x
,
y
)
]
γ
+
k
2
[
I
n
(
x
,
y
)
]
2
γ
+
k
3
[
I
n
(
x
,
y
)
]
3
γ
I_n^c\left( {x,y} \right) = {k_0} + {k_1}{\left[ {{I_n}\left( {x,y} \right)} \right]^\gamma } + {k_2}{\left[ {{I_n}\left( {x,y} \right)} \right]^{2\gamma }} + {k_3}{\left[ {{I_n}\left( {x,y} \right)} \right]^{3\gamma }}
Inc(x,y)=k0+k1[In(x,y)]γ+k2[In(x,y)]2γ+k3[In(x,y)]3γ
由此可见,相机拍摄的光栅图案,存在着高阶(
3
γ
3 \gamma
3γ)的非线性,而光栅图像的非线性表现,表现为图像存在高次谐波,也就是我们常见的水波纹。
03 误差分析
我们实际拍摄到的光栅图案 [48,101,102]:
I
n
c
(
x
,
y
)
=
f
(
I
n
(
x
,
y
)
)
=
a
0
+
∑
k
=
1
K
a
k
cos
(
k
[
φ
(
x
,
y
)
+
δ
n
]
)
I_n^c\left( {x,y} \right) = f\left( {{I_n}\left( {x,y} \right)} \right) = {a_0} + \sum\limits_{k = 1}^K {{a_k}\cos \left( {k\left[ {\varphi \left( {x,y} \right) + {\delta _n}} \right]} \right)}
Inc(x,y)=f(In(x,y))=a0+k=1∑Kakcos(k[φ(x,y)+δn])
其中:
- K K K:最大谐波次数
- a k a_k ak: k k k 次谐波的系数
对这个实际光栅进行解码得到的相位值:
φ
′
(
x
,
y
)
=
−
arctan
[
∑
n
=
0
N
−
1
{
a
0
+
∑
k
=
1
K
a
k
cos
(
k
[
φ
(
x
,
y
)
+
δ
n
]
)
}
sin
(
δ
i
)
∑
n
=
0
N
−
1
{
a
0
+
∑
k
=
1
K
a
k
cos
(
k
[
φ
(
x
,
y
)
+
δ
n
]
)
}
cos
(
δ
i
)
]
\varphi '\left( {x,y} \right) = - \arctan \left[ {\frac{{\sum\limits_{n = 0}^{N - 1} {\left\{ {{a_0} + \sum\limits_{k = 1}^K {{a_k}\cos \left( {k\left[ {\varphi \left( {x,y} \right) + {\delta _n}} \right]} \right)} } \right\}\sin \left( {{\delta _i}} \right)} }}{{\sum\limits_{n = 0}^{N - 1} {\left\{ {{a_0} + \sum\limits_{k = 1}^K {{a_k}\cos \left( {k\left[ {\varphi \left( {x,y} \right) + {\delta _n}} \right]} \right)} } \right\}\cos \left( {{\delta _i}} \right)} }}} \right]
φ′(x,y)=−arctan⎣⎢⎢⎢⎡n=0∑N−1{a0+k=1∑Kakcos(k[φ(x,y)+δn])}cos(δi)n=0∑N−1{a0+k=1∑Kakcos(k[φ(x,y)+δn])}sin(δi)⎦⎥⎥⎥⎤
通常,五次以上的高次谐波对相位计算影响非常小,因此,后续计算中
K
K
K 的取值:
1
−
5
1-5
1−5。而其是由理想相位值
φ
\varphi
φ 与相位误差
Δ
φ
\Delta \varphi
Δφ 构成的,即:
φ
′
(
x
,
y
)
=
φ
(
x
,
y
)
+
Δ
φ
(
x
,
y
)
\varphi '\left( {x,y} \right) = \varphi \left( {x,y} \right) + \Delta \varphi \left( {x,y} \right)
φ′(x,y)=φ(x,y)+Δφ(x,y)
对于标准的三步相移,解码得到的相位值:
φ
′
(
x
,
y
)
=
arctan
[
a
1
sin
[
φ
(
x
,
y
)
]
−
a
2
sin
[
2
φ
(
x
,
y
)
]
+
a
4
sin
[
4
φ
(
x
,
y
)
]
−
a
5
sin
[
5
φ
(
x
,
y
)
]
a
1
cos
[
φ
(
x
,
y
)
]
+
a
2
cos
[
2
φ
(
x
,
y
)
]
−
a
4
cos
[
4
φ
(
x
,
y
)
]
+
a
5
cos
[
5
φ
(
x
,
y
)
]
]
\varphi '\left( {x,y} \right) = \arctan \left[ {\frac{{{a_1}\sin \left[ {\varphi \left( {x,y} \right)} \right] - {a_2}\sin \left[ {2\varphi \left( {x,y} \right)} \right] + {a_4}\sin \left[ {4\varphi \left( {x,y} \right)} \right] - {a_5}\sin \left[ {5\varphi \left( {x,y} \right)} \right]}}{{{a_1}\cos \left[ {\varphi \left( {x,y} \right)} \right] + {a_2}\cos \left[ {2\varphi \left( {x,y} \right)} \right] - {a_4}\cos \left[ {4\varphi \left( {x,y} \right)} \right] + {a_5}\cos \left[ {5\varphi \left( {x,y} \right)} \right]}}} \right]
φ′(x,y)=arctan[a1cos[φ(x,y)]+a2cos[2φ(x,y)]−a4cos[4φ(x,y)]+a5cos[5φ(x,y)]a1sin[φ(x,y)]−a2sin[2φ(x,y)]+a4sin[4φ(x,y)]−a5sin[5φ(x,y)]]
相位误差:
其中:
c
1
c_1
c1 和
c
2
c_2
c2 为常数,且通常
c
1
c_1
c1 要远大于
c
2
c_2
c2,因此相位误差可近似。
其他步长的相移误差同理可得:
由此可推得:当拍摄的光栅图像存在高次谐波时,由:3步相移算法计算的相位,会出现三倍频的周期误差,4步,四倍频。