结构光三维重建这个坑,我断断续续折腾了大概两周。最开始的感受是:论文里格雷码讲得头头是道,相移公式也推得明明白白,但自己动手写代码后才发现,真正把这两个东西配合起来工作,细节才是主角。格雷码负责用一组黑白条纹确定"当前像素对应投影仪的哪根条纹",相移负责在单根条纹内部把相位精确到亚像素,两者一粗一细,刚好互补。这篇就用Python把格雷码加相移这套结构光三维重建的核心链路完整串一遍,把我踩过的坑、改过的方案、以及实测中真正有效的代码片段都摊开讲。
如果你正打算入门结构光三维重建,或者已经看过不少资料但不知道代码从哪里下手,这篇应该对你有用。我会按照实际重建流程组织内容:编码原理、设备准备、格雷码解码、相移计算、相位融合、三角测量,最后是实测避坑。每个环节我都会给出能直接运行的Python核心代码,并解释为什么这么写、为什么这一步非这么做不可。
1. 出发点:为什么非要把格雷码和相移绑在一起用
1.1 先梳理一遍结构光重建的完整链路
结构光三维重建的逻辑可以压缩成一句话:投影仪向物体表面投射一组已知编码图案,相机拍下被物体表面调制后的变形图像,然后根据图像上的灰度变化反推出空间位置关系。想得到三维点云,至少需要三样东西:相机图像上每个像素的坐标、投影仪图像上对应像素的坐标、以及相机和投影仪之间的空间位置关系。
前两样合起来叫对应点匹配,是整个结构光重建的核心。格雷码加相移这套组合,就是专门用来做高精度对应点匹配的。格雷码先给出大范围的绝对索引——物体表面第N根条纹对应投影仪第N个周期;相移再在这个周期内部做精细定位——精确到周期内的相位小数。两者组合之后,每个相机像素都能唯一、连续地对应到投影仪图案上的一个亚像素坐标。
不过这中间有一个经常被忽略的前提:相机拍摄到的图案必须和预处理时生成的图案一一对应。顺序一乱,后面所有编码全部作废。所以我推荐先把主链路印在脑子里,再动手写代码:
- 生成格雷码编码图案和相移图案;
- 按固定顺序投影并同步采集;
- 对格雷码图案做二值化解码,得到条纹索引;
- 对相移图案计算包裹相位;
- 用格雷码索引给每个相移周期定级,把包裹相位展开成绝对相位;
- 用绝对相位建立相机像素与投影仪像素的对应;
- 结合标定参数三角测量,生成三维点云。
整套流程里,融合那一步最容易翻车,后面我会专门讲。
1.2 格雷码的粗定位能力和它的天花板
格雷码的编码逻辑,是把条纹按二进制方式逐位拆开。投影N幅黑白条纹图,每幅图代表一个二进制位,2^N个条纹区域就能用N位二进制码唯一标识。举个例子:投影6幅格雷码图案,就能区分64根条纹区域,条纹越窄,空间分辨率越高。
格雷码和普通二进制编码最大的区别,在于相邻两个码字之间只相差一个比特。这样设计的原因很直接:在条纹边界处,亮度判断因为噪声或物体表面反射问题容易出错,如果相邻码字只差一位,出错时最多落到相邻的条纹区域,不会跳到十万八千里外,后面相移展开时就不会出现大面积的相位跳变。
但是格雷码本身有分辨率天花板。它只能做到"像素级"的条纹定位,想要在单个条纹内部继续细分,格雷码就没辙了。所以必须引入相移。
1.3 相移的高精度和绝对相位缺失问题
相移的思路是投射强度随位置呈正弦变化的条纹图案。通过多幅正弦条纹的相位偏移,可以精确计算出每个像素在正弦周期内部的相位值。四步相移是最经典的形式,用四幅相位依次相差90度的正弦条纹,求一个反正切,就能得到亚像素级别的相位。它在单周期内的精度非常高,噪声影响也小。
但相移有个致命弱点:反正切算出来的相位被压缩在(-π, π]区间内。也就是说,几个相邻的正弦周期在计算结果上是"长得一样"的,无法区分自己位于第几个周期。这就是相位包裹问题。如果不解决周期歧义,最后重建出来的物体表面会出现大量错位条纹,根本无法形成连续曲面。
于是格雷码的价值就体现出来了:格雷码负责告诉每个像素"你现在处于第几个正弦周期",相移负责精确算出"在这个周期里你处于什么亚像素位置"。两者分工不同,缺一不可。
2. 设备、标定与同步,先给重建搭个台子
2.1 最简设备清单和标定常识
结构光重建需要的最低配置并不高:一台工业相机、一台支持外部触发或至少能手动控制的投影仪、一台电脑。相机最好是灰度工业相机,因为格雷码和相移本质上都是对亮度编码,灰度图能避免彩色通道间的串扰。投影仪要求支持纯图案输出,分辨率最好不低于1024像素宽,否则条纹细分不够。
这里必须说清楚:相机和投影仪之间的位姿关系不是自己猜的,需要提前标定。标定思路一般是先固定相机,用棋盘格标定相机的内参和畸变,再通过投影棋盘格或其他编码图案来确定投影仪相对于相机的空间关系。具体公式是相机三维点投影到像素坐标的针孔模型,以及投影仪被当作一个"反向相机"时的投影矩阵。
2.2 相机和投影仪的触发同步,容易翻车的一环
投影仪每切换一幅图案,相机就要拍一张,两边必须严格同步。如果靠手动一帧一帧拍,基本没法用,因为投影仪图像刷新和相机曝光之间存在延迟,稍微错一帧,条纹位置就对不上,解码结果会和实际条纹错位好几个周期。
最稳的方案是硬件同步:电脑把当前图案帧输出到投影仪,同时通过外触发信号让相机开始曝光。如果你的相机不支持外触发,那至少要保证曝光时间远大于投影仪刷新周期,让每帧图像都包含完整的当前条纹图案,然后通过GPIO或串口来控制切换顺序。实测中,我还遇到过投影仪内置的显示缓存导致画面滞后一帧的问题,解决办法是在生成图案序列时,开头插入两幅"空场"图,丢弃不采集,等投影仪刷新稳定后再开始正式采图。
2.3 投影序列设计:一共要投多少张图
以编码64个条纹区域、相移周期为每个区域宽度32像素为例,假设投影仪分辨率是1024像素宽,那么条纹宽度是16像素,需要6位格雷码图案,加上4步相移,总共10幅图。如果物体表面的反光复杂,或追求更高鲁棒性,可以再加上互补格雷码,又要多投好几幅。设计投影序列时,建议把格雷码图案放在前面、相移图案放在后面,顺序固定写进脚本,采集和解码共用同一套顺序。
3. 格雷码解码:把条纹序号从图像里抠出来
3.1 二值化不能只用固定阈值
格雷码解码的第一步,是把图案图像转成0/1序列。很多初学者直接拿一个固定阈值,比如灰度大于128就算亮,然后就会在反光区域或暗角区域翻车。因为结构光投影到物体表面后,不同材质、不同角度的反射率差别很大,同一幅图案里亮区灰度可能从60到240,暗区也可能从10到80,全局阈值很难同时适配整张图。
更靠谱的办法是使用局部自适应阈值。我是用了一种很朴素但有效的做法:对每组格雷码图案,先生成它们的均值图,再以均值图的一定比例作为阈值。具体来说:
import numpy as np def adaptive_threshold(gray_images, k=0.5): # gray_images: 该组所有格雷码图案组成的列表,每张都是二维 numpy 数组 means = np.mean(gray_images, axis=0).astype(np.float32) threshold = means * k binary_maps = [] for img in gray_images: binary = (img > threshold).astype(np.uint8) binary_maps.append(binary) return binary_maps每个像素的阈值来自它在所有格雷码图案中的平均亮度,能够在暗反射区域自动降低阈值,在高光区域自动抬高阈值。实测下来比固定阈值稳定很多。
3.2 格雷码到二进制码的转换逻辑
先解释一下投影端的格雷码序列如何生成。对第 i 根条纹,先用格雷码公式得到它的格雷码:
def gray_encode(n): return n ^ (n >> 1)然后对每一位bit,把对应bit为1的条纹区域置为白色,为0的置为黑色,生成一幅图案。这样拍N幅图,每个像素就能按位组合出一个N位格雷码。接下来要把格雷码还原成二进制码,标准递推关系是:最高位保持不变,之后每一位都是当前格雷码位与上一位二进制位做异或。写成代码:
def gray_to_binary(gray_code_sequence): # gray_code_sequence: 形状为 (N, H, W) 的 0/1 数组 n, h, w = gray_code_sequence.shape binary = np.zeros((h, w), dtype=np.uint32) prev = np.zeros((h, w), dtype=np.uint32) for i in range(n): gray_bit = gray_code_sequence[i].astype(np.uint32) prev = prev ^ gray_bit binary = (binary << 1) | prev return binary解码完成后,binary数组里的每个值就是该像素对应的条纹索引。这一步其实相当快,主要耗时在I/O和图像预处理上。
3.3 解码完成后别急着往下走
格雷码解码结果有个典型问题:条纹边界处会出现大量噪点和毛刺。由于投影图案的亮暗边界在像素尺度上是跳变的,相机拍出来后边缘会有半灰度过渡带,二值化时就会不稳定。边界像素被分配到错误的条纹索引,后面相移展开时就会出现一整行或一整列的周期跳变。
一个缓解手段是对解码结果做中值滤波,但要注意:中值滤波的窗口不要太大,否则会把细窄条纹的细节抹掉。我用的是3x3窗口的中值滤波,只在明显孤立的毛刺区域做修正,不破坏条纹边界本身。更彻底的办法会在后面融合那节单独讨论。
4. 相移计算:用四步公式抠出亚像素相位
4.1 四步相移的数学原理
四步相移法的原理不复杂:四幅正弦条纹图案,初始相位依次相差90度。第k幅图案的灰度可以写成:
I1 = A + B cos(φ) I2 = A + B cos(φ + 90°) I3 = A + B cos(φ + 180°) I4 = A + B cos(φ + 270°)
其中A是平均亮度,B是条纹对比度,φ就是我们要的相位。对四幅图做数学变形,可以消掉A和B,最终得到:
φ = arctan((I4 - I2) / (I1 - I3))
注意这里的除法直接用会丢失象限信息,所以一定要用atan2函数,它不仅能保证相位落在(-π, π]区间内,还能正确处理分子分母为负的情况。Python代码很简短:
def compute_wrapped_phase(I1, I2, I3, I4): numerator = I4.astype(np.float32) - I2.astype(np.float32) denominator = I1.astype(np.float32) - I3.astype(np.float32) phase = np.arctan2(numerator, denominator) return phase4.2 包裹相位到底是什么
上面算出来的phase,是一个位于(-π, π]之间的值。如果一个物体表面的深度变化跨过了多个正弦周期,相位值就会在超过π的时候突然跳到-π,形成锯齿状的"包裹相位"图。直观理解就是:你明明知道每根条纹对应的连续相位应该是0到多个2π递增,但计算出来却总是在π附近反复折叠。
这个折叠不是噪声,而是反正切函数的数学特性。所以后续一定要做相位展开。相位展开的思路是:对每个像素,先确定它处于第几个周期(格雷码已经算出来了),然后把这个周期序号乘上2π,再加回包裹相位,就得到了连续的绝对相位。
4.3 为什么要四步而不是三步
三步相移相比四步也能算相位,只需要三幅图,公式是:
φ = arctan(√3 * (I2 - I3) / (2*I1 - I2 - I3))
优点是少投一幅图,采集时间更短;缺点是对平均分量和伽马畸变更敏感。四步相移公式中分母自带一些抵消效果,抗噪声能力更好,重建出来的条纹平滑度也更高。我实测下来,动态场景用三步,静态精细扫描用四步。既然这篇聊的是能把精度做高的一整套方案,我就以四步为例。
5. 融合与展开:相位和格雷码如何配合
5.1 格雷码条纹边界导致的相位跳变问题
假设投影仪一行有1024个像素,相移周期宽度是32像素,那么一行里共有32个周期。如果格雷码完全准确,每个像素的周期索引范围是0到31,相位展开公式就是:
absolute_phase = wrapped_phase + 2π * period_index
但格雷码解码并不是完美的。尤其在条纹边界的像素上,解码出的索引可能和相邻像素差一档,于是展开后的绝对相位连续表面会硬生生裂成两半。这个问题在重建物体表面时表现为明显的台阶或横线。
5.2 互补格雷码:磨平边界误差的常见方案
一种非常实用的手段是"互补格雷码"。思路是:投影两组格雷码图案,第二组比第一组移动半个条纹宽度。这样一来,第一组的条纹边界和第二组的条纹边界错开了半个条纹位置。对每个像素分别用两组格雷码解码,如果两个索引一致,就取该值;如果不一致,说明该像素落在某个条纹边界附近,直接标记为无效或保守处理。
实测下来,互补格雷码能消灭绝大多数边界跳变,代价是投影图案数量几乎翻倍。如果对采集速度要求高,也可以用第一组解码加上相位连续性约束来做后处理,但那属于修bug的思路,远不如一开始就把边界处理掉省心。
5.3 融合代码:从周期索引到绝对相位
合并格雷码索引和包裹相位的核心代码其实很简洁:
def unwrap_phase_with_gray(wrapped_phase, period_index, period_per_image): # period_index 是格雷码解码出的周期序号,从 0 开始 # wrapped_phase 在 [-pi, pi) 之间 # 投影时周期宽度可能是 32 像素,格雷码条纹宽度也是 32 像素,所以索引直接对应 absolute_phase = wrapped_phase + 2.0 * np.pi * period_index.astype(np.float32) return absolute_phase这里有一个隐藏前提:格雷码条纹宽度必须和相移周期宽度一致。如果格雷码一根条纹覆盖64像素,而相移周期只有32像素,那就得先把周期索引乘2再融合,否则会发生错位。很多新手在这里踩坑,根源就是这个宽度比例没对齐。
融合后得到的absolute_phase沿着条纹方向单调递增,可以作为相机图像到投影仪坐标的桥梁。如果物体表面深度变化不大,不需要显示地做二维相位展开,只要保证每行/每列都能正确解码和融合就行。
6. 三角测量:从对应点到三维坐标
6.1 把投影仪当作一台逆相机
三角测量的输入,是相机图像坐标(u_c, v_c)和它对应的投影仪坐标(u_p, v_p)。相机坐标可以直接用图像像素位置表示;投影仪坐标怎么来?绝对相位就是投影仪坐标的编码。如果投影的是竖直条纹,那么相位值和投影仪像素列坐标是一一映射关系:
u_p = (absolute_phase / (2π)) * (W_p / num_periods)
其中W_p是投影仪分辨率宽度,num_periods是整幅投影图案上的周期总数。竖直条纹只对列编码,所以v_p一般不需要靠相移求,可以直接取0到高度之间的一系列值,或者用另一组水平条纹单独编码。
6.2 标定矩阵的三角化方程
标定之后,相机和投影仪分别对应一个3x4投影矩阵P_c和P_p。对空间中的三维点X,它在相机和投影仪上的齐次坐标分别满足:
λ_c * (u_c, v_c, 1)^T = P_c * (X, 1)^T
λ_p * (u_p, v_p, 1)^T = P_p * (X, 1)^T
两个方程一共四个线性约束,未知数只有X的三个分量,属于超定方程。常规做法是用最小二乘解:
def triangulate(u_c, v_c, u_p, v_p, P_c, P_p): # P_c, P_p 都是 3x4 的投影矩阵,numpy 数组 A = np.zeros((4, 4), dtype=np.float32) b = np.zeros((4,), dtype=np.float32) A[0, :3] = u_c * P_c[2, :3] - P_c[0, :3] A[1, :3] = v_c * P_c[2, :3] - P_c[1, :3] A[2, :3] = u_p * P_p[2, :3] - P_p[0, :3] A[3, :3] = v_p * P_p[2, :3] - P_p[1, :3] b[0] = P_c[0, 3] - u_c * P_c[2, 3] b[1] = P_c[1, 3] - v_c * P_c[2, 3] b[2] = P_p[0, 3] - u_p * P_p[2, 3] b[3] = P_p[1, 3] - v_p * P_p[2, 3] X, res, rank, sv = np.linalg.lstsq(A, b, rcond=None) return X这段代码看起来短,但每一行都是从投影矩阵方程推出来的,不是硬背出来的。理解它的关键在于:每一条像素对应关系,其实是一条三维射线,两个射线的交点就是三维点。由于噪声存在,两条射线不一定严格相交,最小二乘会求出一个误差最小的近似交点。
6.3 点云生成与后处理
把相机图像里的每一个有效像素都做一次三角测量,就得到一张点云。不过这么做数据量不小,而且很多像素是无效的,比如没有投影覆盖的区域或反射过强的区域。实际处理时我会先做一次掩膜,只对格雷码解码结果稳定、且相移图亮度在合理范围内的像素进行三角测量,能省掉一大半无效计算。
点云生成后还要做一步去噪。简单做法是统计每个点的邻域距离,距离明显大于周围平均值的点直接删除,再对点云法线做一次一致性修正。结构光的误差分布和激光雷达不太一样,主要在条纹边界处会出现少量离群点,所以去噪阈值不用设太小,重点还是在"条纹边界处理"这一关。
7. 实测避坑与调优笔记
7.1 解码条纹错位怎么办
我遇到过最典型的错位是:重建出的物体表面每隔固定间隔出现一道台阶。检查后发现不是相移问题,而是格雷码解码时有个别条纹索引完全偏了一整级。原因很有意思:我生成的格雷码图案顺序和解码时读取的顺序不一致,导致每个像素都被平移了一个位。后来我在代码里加了校验:对物体表面点亮一点的平面区域,绝对相位沿水平方向应该单调递增,如果出现逆转,说明解码顺序反了或者条纹宽度和周期宽度没对齐。
7.2 反光表面和高光区域的处理
高光区域是格雷码和相移共同的大敌。物体表面某个点反射过强,相移图案的灰度会直接饱和,算出的相位完全失真。我实测的解决办法有三种:一是投影前先拍一张环境光图,在正式计算时减去环境光;二是曝光时间尽量压低,以暗部不丢失为底线,避免亮部饱和;三是在掩膜阶段把高光像素直接剔除,因为它们恢复出来的深度本来就不可信。
还有一种情况是物体表面存在互反射,比如金属凹槽或相邻平面,导致投影条纹被二次反射到不该出现的地方。这类问题单纯靠调整格雷码和相移已经很难解决,需要借助深度学习或更复杂的光学方案,不在这篇讨论范围内。
7.3 性能调优:从跑几分钟到几秒钟
纯Python跑完整套流程,速度确实不够理想,瓶颈主要在逐像素循环上。我的调优顺序是:
- 所有图像处理先用numpy向量化,避免Python循环遍历像素;
- 格雷码解码的位操作全部用numpy位运算完成;
- 三角测量不要对整个图像做遍历循环,而是用矩阵形式一次性处理所有坐标点;
- 如果内存允许,可以把多幅图像堆叠成三维数组,利用并行批处理。
优化完之后,1024x768分辨率的点云生成,在一次静态模式下大约两秒左右能跑完,足够做离线扫描和逆向建模。如果还想更快,就得考虑GPU或CUDA版的三角测量了,但那已经有点超出"手撕代码"的范畴。
折腾完这一整轮,我的感受是:格雷码和相移这套组合,难的不是看懂公式,而是把公式、编码顺序、设备同步、标定参数全部揉到一起不打架。只要每一步都按照"粗定位到精定位、周次到亚像素相位、对应关系到三维坐标"这条逻辑走下来,结构光三维重建其实没有想象中那么神秘。最后再提醒一句:遇到任何奇怪的条纹重影,先查编码顺序,先查条纹宽度和周期宽度的对齐关系,八成问题都出在那些看起来最不起眼的地方。