简介面向光学测量、机器视觉及数字全息干涉方向的学习者该仿真包提供相移轮廓术的代码实现。内容围绕四步相移算法覆盖条纹图案生成、相位解算、相位解包裹与高度恢复等关键环节帮助理解相位与物体表面高度之间的关系。压缩包共5个文件包括两个m脚本主程序main.m与解包裹函数jiebaoguo.m、一幅tif重构结果图、一份PDF说明文档及一张jpg物体样图整体仅206KB轻量易用。已有992人学习下载适合刚接触结构光三维测量的学生或工程师快速上手。通过运行代码可直观观察条纹投影到三维形貌重建的完整流程同时结合PDF文档理解四步相移算法的数学推导为后续开展精密测量实验或工程应用打下基础。1. 相移轮廓术仿真实验到底在仿什么相移轮廓术Phase Shifting Profilometry仿真的重点不是把条纹画得漂亮而是给自己一个真值已知的三维测量闭环。真机系统里投影仪畸变、相机噪声、Gamma 非线性、环境光混在一起算法偏了很难定位问题出在哪一环仿真里高度场是亲手定义的条纹按公式从这个高度场生成于是包裹相位、解包裹、相位-高度映射每一步都能和真值逐像素对账。压缩包常见形态是 Python 或 MATLAB 脚本输入合成相移条纹输出高度图或三维点云适合刚接触光栅投影三维测量的学生也适合想先把算法流程跑通再上真机验证的工程师。仿真的最终产出不是渲染图而是一张误差分布表。2. 相移轮廓术仿真第一步条纹生成与采样误差建模2.1 仿真系统里的投影-相机模型与 4 个关键参数不需要真的建模投影仪光路常见做法是把它当成逆相机投影仪投出已知相位的正弦光栅相机看到的是被物体高度调制后的同一光栅。把相机像素坐标记为 (x,y)第 k 帧的强度写成I_k(x,y) A(x,y) B(x,y)cos[φ(x,y) 2πk/N]。其中 A 是背景光B 是调制度也就是条纹对比度N 是相移步数。待求相位 φ 包含两项沿 x 单调增加的载频项 2πx/p以及高度引起的相位偏移 2πh/Λ。Λ 是等效波长物理含义就是相位变化 2π 时高度变化多少。仿真里这 4 个参数直接决定后面所有环节的难度先列成表参数符号仿真里的取值影响条纹周期p32 px越小载频越高解包裹越敏感相移步数N4备选 3/6/12噪声与谐波抑制能力背景光/调制度A / B0.5 / 0.5决定灰度动态范围和信噪比等效波长Λ10高度与相位的换算比例我一般建议先固定 p32、N4 跑通全流程再动其它参数。p 选得太小会让相邻像素相位差逼近 π这是空间解包裹的硬上限B 调太低则包裹相位里的噪声会大到解包裹跳线的程度。这两个值是仿真实验里最常见的翻车点。2.2 用 NumPy 生成 N 步相移条纹与物体高度场生成部分用 NumPy 就够了不需要仿真软件。下面这段生成一个 512×512 的半圆帽物体并合成 4 步相移条纹import numpy as np rng np.random.default_rng(2024) # 固定种子保证每次结果可比 h, w 512, 512 y, x np.mgrid[0:h, 0:w].astype(np.float64) # 半圆帽物体直径 240 px高度 0~16与 lam_eq 同单位 r np.sqrt((x - 256.0)**2 (y - 256.0)**2) mask r 120.0 h_true np.zeros((h, w)) h_true[mask] 8.0 * np.cos(np.pi * r[mask] / 240.0) 8.0 pitch 32 # 条纹周期单位像素 lam_eq 10.0 # 等效波长相位 2π 对应高度 10 plane_k 2.0 * np.pi / pitch # 参考平面相位沿 x 的斜率 phi plane_k * x 2.0 * np.pi * h_true / lam_eq N 4 fringes [] for k in range(N): delta 2.0 * np.pi * k / N fringe 0.5 0.5 * np.cos(phi delta) # A0.5, B0.5 fringes.append(fringe) fringes np.stack(fringes) # 形状 (N, h, w)注意两点。一是相位里先写载频项再写高度项载频项让包裹相位在 x 方向有规律地跳 2π正好用来验证解包裹能不能把这些跳变修回来二是 A 和 B 都取 0.5让灰度落在 0~1 且不饱和后面加噪声、量化才有余量。如果你把灰度直接按 0~255 生成噪声尺度也要跟着乘 255这是很常见的取值混乱。2.3 给理想条纹加三类真实感误差理想条纹直接进算法会得到几乎完美的结果看不出算法边界。我一般会分三步污染先对条纹做高斯低通等效投影仪离焦再加传感器加性高斯噪声最后按 8 bit 量化。顺序不要颠倒它对应光学模糊、传感器采样、ADC 量化的物理顺序。from scipy.ndimage import gaussian_filter blur_sigma 1.2 # 离焦半径像素 noise_sigma 0.01 # 噪声标准差0~1 灰度尺度 for k in range(N): fringes[k] gaussian_filter(fringes[k], sigmablur_sigma) fringes[k] rng.normal(0.0, noise_sigma, (h, w)) fringes_q np.clip(np.round(fringes * 255.0), 0, 255).astype(np.uint8)调 blur_sigma 时注意一个边界离焦过重会让调制度 B 下降物体边缘的相位会漂移这种误差是系统性的增大 N 救不回来只能靠边缘剔除或去卷积。量化这一步在仿真里最容易忽略但 8 bit 量化本身会引入约 1/√12 ≈ 0.29 个灰度级的强度量化噪声换算成相位还要除以调制度和 √(N/2)肉眼看不见但误差预算里要留这笔。3. 相移轮廓术核心参数N 步相移公式与包裹相位提取3.1 用 atan2 从 N 步相移条纹提取包裹相位N 步相移的提取公式是所有相移轮廓术算法的地基。对前面生成的 I_k先分别对相移量 δ_k2πk/N 做正弦、余弦加权求和s(x,y)Σ_k I_k sin δ_kc(x,y)Σ_k I_k cos δ_k。包裹相位就是 φ_w -atan2(s, c)取值范围 (-π, π]。这里的负号来自我们采用的相移符号约定不同文献可能把相移写成 φ-δ_k公式里的负号会反过来这是后续所有现象是否上下颠倒的总根源。写成代码def phase_shift_reconstruct(fringes): N 步相移还原返回包裹相位(rad)与调制度质量图。 N fringes.shape[0] k np.arange(N, dtypenp.float64) delta 2.0 * np.pi * k / N s np.sum(fringes * np.sin(delta), axis0) c np.sum(fringes * np.cos(delta), axis0) phi_w -np.arctan2(s, c) # 包裹到 [-pi, pi) B 2.0 / N * np.sqrt(s**2 c**2) # 调制度兼作质量图 return phi_w, B phi_w, B phase_shift_reconstruct(fringes_q) # 直接吃第 2 节量化后的图为什么要用 atan2 而不是 atan因为 atan 的返回值只有两象限s、c 同号相除时会把 (-π, -π/2) 和 (π/2, π) 这两个象限弄混但 atan2 能根据 s、c 的符号还原完整象限直接给出 2π 范围内的角度。B 这一路输出是白拿的副产品它等于条纹调制度物体表面反射率低、离焦严重的区域 B 会明显变小后面解包裹就用它当质量图。3.2 N 取 3、4、6、12 时噪声与谐波误差怎么变N 的选择直接决定结果质量和采集成本。随机噪声方面相位误差标准差近似与 1/√N 成正比步数越多越抗噪。谐波方面有个简单判据N 步相移只能抑制阶次 m 满足 m ≢ ±1 (mod N) 的谐波其余阶次会原样泄漏进相位误差。投影仪 Gamma 非线性主要产生二次谐波所以用 3 步相移时二次谐波正好在泄漏区间误差最明显4 步相移则把偶次谐波全压掉只留奇次。整理成常用对照表步数 N相对相位噪声泄漏的谐波阶次典型场景31.002、4、5、7…高速动态测量40.873、5、7…奇次通用静态测量60.715、7、11…中等噪声、兼顾速度120.5011、13…高精度静态测量我在做仿真对比时固定其它参数只改 N把几次重复实验的相位误差标准差画在同一张图上基本都能看到 1/√N 的下降趋势。如果曲线明显偏离先去查条纹里是不是混了谐波或者噪声不是高斯的——比如量化噪声在 B 极小时不再是白噪声这时加大 N 的收益会递减。3.3 包裹相位阶段常犯的 3 个错误第一个是符号反了。φ_w 差一个负号时重建出来的高度场会整个凹变凸而且由于载频项的存在你甚至不会立刻发现方向反了。验证方法很简单取物体中心一列像素把 φ_w 与输入的 phi_obj 画在一起看趋势是否一致。第二个是拿 uint8 图直接算。亮度饱和或截断的区域s、c 的线性关系被破坏atan2 出来的相位会带上系统误差仿真里建议保留浮点分支做谐波对比实验时单独开关量化这一项。第三个是忽略 B 阈值。物体边缘、阴影区的 B 趋近于零那里的包裹相位纯属噪声解包裹时却会把误差沿等值线传播开。正确做法是把 B 0.3 的像素先打上掩膜不进解包裹最后重建也不统计这些点。4. 相移轮廓术解包裹与相位-高度映射4.1 质量图引导解包裹优先展开调制高的像素包裹相位在 (-π, π] 内跳变空间解包裹的思路是逐像素比较相邻相位遇到大于 π 的跳变就补一个 2π。质量图引导算法是其中最稳的一种从调制度 B 最高的像素出发用优先队列总是先展开邻居已展开、自身质量最高的像素把 2π 跳变逐步从高质量区向低质量区推开。import heapq def quality_guided_unwrap(phi, q): 质量图引导解包裹q 为调制度质量图。 hgt, wdt phi.shape visited np.zeros((hgt, wdt), dtypebool) unwrapped np.zeros_like(phi) si, sj np.unravel_index(np.argmax(q), phi.shape) visited[si, sj] True unwrapped[si, sj] phi[si, sj] heap [] for di, dj in ((1, 0), (-1, 0), (0, 1), (0, -1)): ni, nj si di, sj dj if 0 ni hgt and 0 nj wdt: heapq.heappush(heap, (-q[ni, nj], ni, nj)) while heap: _, i, j heapq.heappop(heap) if visited[i, j]: continue visited[i, j] True vals [] for di, dj in ((1, 0), (-1, 0), (0, 1), (0, -1)): ni, nj i di, j dj if 0 ni hgt and 0 nj wdt and visited[ni, nj]: vals.append(unwrapped[ni, nj]) base float(np.mean(vals)) # 已展开邻域的平均相位 k round((base - phi[i, j]) / (2.0 * np.pi)) # 最近整数个 2π unwrapped[i, j] phi[i, j] 2.0 * np.pi * k for di, dj in ((1, 0), (-1, 0), (0, 1), (0, -1)): ni, nj i di, j dj if 0 ni hgt and 0 nj wdt and not visited[ni, nj]: heapq.heappush(heap, (-q[ni, nj], ni, nj)) return unwrapped优先队列里用负质量当键pop 出来就是质量最高的待展开像素。每展开一个像素用已展开邻域相位均值做基准计算当前包裹相位与基准相差的 2π 整数倍并补上这就是展开的本质。这个实现里没有处理遮挡和断崖真实工程还要加掩膜与残差点标记但作为仿真实验的基准已经够用。跑完把 unwrapped 画成伪彩图物体区域颜色应该连续且平滑如果出现一条条接缝说明有跳点没补上。4.2 相位差换算高度等效波长模型与真值对比解包裹得到的绝对相位里还叠着载频项要还原高度得先去掉它。仿真里最干净的做法是同样跑一遍零高度参考平面两幅展开相位相减Δφ(x,y) φ_obj(x,y) - φ_ref(x,y)h(x,y) Λ·Δφ(x,y) / (2π)。这个线性关系只在等效波长模型成立时准确。真实系统里投影仪和相机光轴不平行Δφ 与 h 是接近双曲的非线性关系标定时要拟合高阶多项式仿真实验里没这个必要但你要清楚线性模型是近似。下面把整条流水线封装成一个函数物体和参考平面各跑一遍def run_sim(height_map, seed0): 输入高度场输出展开相位和调制度质量图。 rng np.random.default_rng(seed) phi plane_k * x 2.0 * np.pi * height_map / lam_eq fr [0.5 0.5 * np.cos(phi 2.0 * np.pi * k / N) for k in range(N)] fr np.stack(fr) for k in range(N): fr[k] gaussian_filter(fr[k], sigmablur_sigma) # 参数来自 2.3 节 fr[k] rng.normal(0.0, noise_sigma, (h, w)) fr np.clip(np.round(fr * 255.0), 0, 255).astype(np.uint8) phi_w, B phase_shift_reconstruct(fr) return quality_guided_unwrap(phi_w, B), B phi_obj, B_obj run_sim(h_true, seed1) phi_ref, B_ref run_sim(np.zeros_like(h_true), seed2) h_est (phi_obj - phi_ref) * lam_eq / (2.0 * np.pi)两个坑在这里最容易出现。一是全局 2π 常数解包裹只保证相对连续物体和参考平面的展开相位可能差一个 2π 的整数倍导致整幅高度偏移一个 Λ。仿真里如果发现 h_est 整体偏高或偏低先检查这个而不是去调滤波参数。二是 run_sim 用到的 plane_k、x、blur_sigma 都来自第 2 节定义的全局环境封装函数时不要另起一套参数否则两次仿真的空间基准就对不齐。提示先做一步零高度参考平面残差检查展开相位相减后残差均值应接近 0若接近 ±2π就是全局常数漂移。4.3 用 RMSE 和边界掩膜量化重建误差仿真实验的最后一环也是它区别于真机调参的地方和真值 h_true 逐像素比大小。valid mask (B_obj 0.3) # 物体内部且调制度足够 err h_est - h_true rmse np.sqrt(np.mean(err[valid] ** 2)) print(f有效像素 {valid.sum()}, RMSE {rmse:.4f})掩膜必须加原因有两层物体外面的像素没有真值可对而物体边缘因离焦模糊B 骤降误差会大一个数量级不剔除的话 RMSE 会被边缘主导看不出算法本身的水平。想看误差的空间分布就把 err 用有效掩膜画热力图通常边缘一圈亮、内部接近零如果内部出现周期性亮带那就是谐波泄漏或解包裹跳线和随机噪声的形状完全不同。各处理阶段的误差来源整理如下阶段主要误差来源仿真里怎么看条纹生成离焦、量化、加性噪声对比加噪前后 B 的直方图包裹相位谐波泄漏、N 过小用 N12 的结果当近似真值解包裹梯度超过 π/px、低质量区检查展开相位的接缝高度映射全局 2π 常数漂移参考平面区域残差均值5. 相移轮廓术仿真的收尾技巧固定种子、中间结果与 git 管理5.1 固定随机种子并逐级落盘仿真实验改一个参数就要重跑一次如果随机种子每次不一样两轮结果的差异分不清是参数引起的还是噪声引起的。我习惯在 run_sim 里把 seed 作为显式参数传入做对比时固定 seed 只改被测参数同时把包裹相位、展开相位、高度估计逐级用 np.save 落盘排查时不用把整条链路重跑一遍。建议的文件布局是 data/ 放条纹与相位中间结果out/ 放高度图和误差图脚本根目录只留源码。5.2 解压校验、EOCD 报错与 git 管理拿到压缩包先做完整性校验再动手改代码。unzip -t 会逐文件测试 CRC如果报 could not find EOCD说明文件没下完整直接用下载工具重试别浪费时间找所谓修复工具。遇到 .z01 分卷或报 CRC 错误的文件用 7z x 统一处理更省事。unzip -t 相移轮廓术仿真实验.zip unzip 相移轮廓术仿真实验.zip -d psp_sim cd psp_sim git init git add -A git commit -m psp simulation baseline解压后第一件事就是 git init。仿真的参数一改结果就变没有版本记录的话一两周后你根本说不清哪张误差表对应哪套参数。每次调完参提交一次commit message 里写清楚改了 p、N 还是噪声对比实验就能随时回退。最后一个验证技巧把展开相位减去载频项 plane_k·x得到的残差应该是平滑的物体相位如果残差图上有锯齿状条纹说明某个像素的 2π 补跳补错了优先检查 B 最低的那片区域。本文还有配套的精品资源点击获取
相移轮廓术仿真实验：条纹生成、相位解包裹与误差量化
简介面向光学测量、机器视觉及数字全息干涉方向的学习者该仿真包提供相移轮廓术的代码实现。内容围绕四步相移算法覆盖条纹图案生成、相位解算、相位解包裹与高度恢复等关键环节帮助理解相位与物体表面高度之间的关系。压缩包共5个文件包括两个m脚本主程序main.m与解包裹函数jiebaoguo.m、一幅tif重构结果图、一份PDF说明文档及一张jpg物体样图整体仅206KB轻量易用。已有992人学习下载适合刚接触结构光三维测量的学生或工程师快速上手。通过运行代码可直观观察条纹投影到三维形貌重建的完整流程同时结合PDF文档理解四步相移算法的数学推导为后续开展精密测量实验或工程应用打下基础。1. 相移轮廓术仿真实验到底在仿什么相移轮廓术Phase Shifting Profilometry仿真的重点不是把条纹画得漂亮而是给自己一个真值已知的三维测量闭环。真机系统里投影仪畸变、相机噪声、Gamma 非线性、环境光混在一起算法偏了很难定位问题出在哪一环仿真里高度场是亲手定义的条纹按公式从这个高度场生成于是包裹相位、解包裹、相位-高度映射每一步都能和真值逐像素对账。压缩包常见形态是 Python 或 MATLAB 脚本输入合成相移条纹输出高度图或三维点云适合刚接触光栅投影三维测量的学生也适合想先把算法流程跑通再上真机验证的工程师。仿真的最终产出不是渲染图而是一张误差分布表。2. 相移轮廓术仿真第一步条纹生成与采样误差建模2.1 仿真系统里的投影-相机模型与 4 个关键参数不需要真的建模投影仪光路常见做法是把它当成逆相机投影仪投出已知相位的正弦光栅相机看到的是被物体高度调制后的同一光栅。把相机像素坐标记为 (x,y)第 k 帧的强度写成I_k(x,y) A(x,y) B(x,y)cos[φ(x,y) 2πk/N]。其中 A 是背景光B 是调制度也就是条纹对比度N 是相移步数。待求相位 φ 包含两项沿 x 单调增加的载频项 2πx/p以及高度引起的相位偏移 2πh/Λ。Λ 是等效波长物理含义就是相位变化 2π 时高度变化多少。仿真里这 4 个参数直接决定后面所有环节的难度先列成表参数符号仿真里的取值影响条纹周期p32 px越小载频越高解包裹越敏感相移步数N4备选 3/6/12噪声与谐波抑制能力背景光/调制度A / B0.5 / 0.5决定灰度动态范围和信噪比等效波长Λ10高度与相位的换算比例我一般建议先固定 p32、N4 跑通全流程再动其它参数。p 选得太小会让相邻像素相位差逼近 π这是空间解包裹的硬上限B 调太低则包裹相位里的噪声会大到解包裹跳线的程度。这两个值是仿真实验里最常见的翻车点。2.2 用 NumPy 生成 N 步相移条纹与物体高度场生成部分用 NumPy 就够了不需要仿真软件。下面这段生成一个 512×512 的半圆帽物体并合成 4 步相移条纹import numpy as np rng np.random.default_rng(2024) # 固定种子保证每次结果可比 h, w 512, 512 y, x np.mgrid[0:h, 0:w].astype(np.float64) # 半圆帽物体直径 240 px高度 0~16与 lam_eq 同单位 r np.sqrt((x - 256.0)**2 (y - 256.0)**2) mask r 120.0 h_true np.zeros((h, w)) h_true[mask] 8.0 * np.cos(np.pi * r[mask] / 240.0) 8.0 pitch 32 # 条纹周期单位像素 lam_eq 10.0 # 等效波长相位 2π 对应高度 10 plane_k 2.0 * np.pi / pitch # 参考平面相位沿 x 的斜率 phi plane_k * x 2.0 * np.pi * h_true / lam_eq N 4 fringes [] for k in range(N): delta 2.0 * np.pi * k / N fringe 0.5 0.5 * np.cos(phi delta) # A0.5, B0.5 fringes.append(fringe) fringes np.stack(fringes) # 形状 (N, h, w)注意两点。一是相位里先写载频项再写高度项载频项让包裹相位在 x 方向有规律地跳 2π正好用来验证解包裹能不能把这些跳变修回来二是 A 和 B 都取 0.5让灰度落在 0~1 且不饱和后面加噪声、量化才有余量。如果你把灰度直接按 0~255 生成噪声尺度也要跟着乘 255这是很常见的取值混乱。2.3 给理想条纹加三类真实感误差理想条纹直接进算法会得到几乎完美的结果看不出算法边界。我一般会分三步污染先对条纹做高斯低通等效投影仪离焦再加传感器加性高斯噪声最后按 8 bit 量化。顺序不要颠倒它对应光学模糊、传感器采样、ADC 量化的物理顺序。from scipy.ndimage import gaussian_filter blur_sigma 1.2 # 离焦半径像素 noise_sigma 0.01 # 噪声标准差0~1 灰度尺度 for k in range(N): fringes[k] gaussian_filter(fringes[k], sigmablur_sigma) fringes[k] rng.normal(0.0, noise_sigma, (h, w)) fringes_q np.clip(np.round(fringes * 255.0), 0, 255).astype(np.uint8)调 blur_sigma 时注意一个边界离焦过重会让调制度 B 下降物体边缘的相位会漂移这种误差是系统性的增大 N 救不回来只能靠边缘剔除或去卷积。量化这一步在仿真里最容易忽略但 8 bit 量化本身会引入约 1/√12 ≈ 0.29 个灰度级的强度量化噪声换算成相位还要除以调制度和 √(N/2)肉眼看不见但误差预算里要留这笔。3. 相移轮廓术核心参数N 步相移公式与包裹相位提取3.1 用 atan2 从 N 步相移条纹提取包裹相位N 步相移的提取公式是所有相移轮廓术算法的地基。对前面生成的 I_k先分别对相移量 δ_k2πk/N 做正弦、余弦加权求和s(x,y)Σ_k I_k sin δ_kc(x,y)Σ_k I_k cos δ_k。包裹相位就是 φ_w -atan2(s, c)取值范围 (-π, π]。这里的负号来自我们采用的相移符号约定不同文献可能把相移写成 φ-δ_k公式里的负号会反过来这是后续所有现象是否上下颠倒的总根源。写成代码def phase_shift_reconstruct(fringes): N 步相移还原返回包裹相位(rad)与调制度质量图。 N fringes.shape[0] k np.arange(N, dtypenp.float64) delta 2.0 * np.pi * k / N s np.sum(fringes * np.sin(delta), axis0) c np.sum(fringes * np.cos(delta), axis0) phi_w -np.arctan2(s, c) # 包裹到 [-pi, pi) B 2.0 / N * np.sqrt(s**2 c**2) # 调制度兼作质量图 return phi_w, B phi_w, B phase_shift_reconstruct(fringes_q) # 直接吃第 2 节量化后的图为什么要用 atan2 而不是 atan因为 atan 的返回值只有两象限s、c 同号相除时会把 (-π, -π/2) 和 (π/2, π) 这两个象限弄混但 atan2 能根据 s、c 的符号还原完整象限直接给出 2π 范围内的角度。B 这一路输出是白拿的副产品它等于条纹调制度物体表面反射率低、离焦严重的区域 B 会明显变小后面解包裹就用它当质量图。3.2 N 取 3、4、6、12 时噪声与谐波误差怎么变N 的选择直接决定结果质量和采集成本。随机噪声方面相位误差标准差近似与 1/√N 成正比步数越多越抗噪。谐波方面有个简单判据N 步相移只能抑制阶次 m 满足 m ≢ ±1 (mod N) 的谐波其余阶次会原样泄漏进相位误差。投影仪 Gamma 非线性主要产生二次谐波所以用 3 步相移时二次谐波正好在泄漏区间误差最明显4 步相移则把偶次谐波全压掉只留奇次。整理成常用对照表步数 N相对相位噪声泄漏的谐波阶次典型场景31.002、4、5、7…高速动态测量40.873、5、7…奇次通用静态测量60.715、7、11…中等噪声、兼顾速度120.5011、13…高精度静态测量我在做仿真对比时固定其它参数只改 N把几次重复实验的相位误差标准差画在同一张图上基本都能看到 1/√N 的下降趋势。如果曲线明显偏离先去查条纹里是不是混了谐波或者噪声不是高斯的——比如量化噪声在 B 极小时不再是白噪声这时加大 N 的收益会递减。3.3 包裹相位阶段常犯的 3 个错误第一个是符号反了。φ_w 差一个负号时重建出来的高度场会整个凹变凸而且由于载频项的存在你甚至不会立刻发现方向反了。验证方法很简单取物体中心一列像素把 φ_w 与输入的 phi_obj 画在一起看趋势是否一致。第二个是拿 uint8 图直接算。亮度饱和或截断的区域s、c 的线性关系被破坏atan2 出来的相位会带上系统误差仿真里建议保留浮点分支做谐波对比实验时单独开关量化这一项。第三个是忽略 B 阈值。物体边缘、阴影区的 B 趋近于零那里的包裹相位纯属噪声解包裹时却会把误差沿等值线传播开。正确做法是把 B 0.3 的像素先打上掩膜不进解包裹最后重建也不统计这些点。4. 相移轮廓术解包裹与相位-高度映射4.1 质量图引导解包裹优先展开调制高的像素包裹相位在 (-π, π] 内跳变空间解包裹的思路是逐像素比较相邻相位遇到大于 π 的跳变就补一个 2π。质量图引导算法是其中最稳的一种从调制度 B 最高的像素出发用优先队列总是先展开邻居已展开、自身质量最高的像素把 2π 跳变逐步从高质量区向低质量区推开。import heapq def quality_guided_unwrap(phi, q): 质量图引导解包裹q 为调制度质量图。 hgt, wdt phi.shape visited np.zeros((hgt, wdt), dtypebool) unwrapped np.zeros_like(phi) si, sj np.unravel_index(np.argmax(q), phi.shape) visited[si, sj] True unwrapped[si, sj] phi[si, sj] heap [] for di, dj in ((1, 0), (-1, 0), (0, 1), (0, -1)): ni, nj si di, sj dj if 0 ni hgt and 0 nj wdt: heapq.heappush(heap, (-q[ni, nj], ni, nj)) while heap: _, i, j heapq.heappop(heap) if visited[i, j]: continue visited[i, j] True vals [] for di, dj in ((1, 0), (-1, 0), (0, 1), (0, -1)): ni, nj i di, j dj if 0 ni hgt and 0 nj wdt and visited[ni, nj]: vals.append(unwrapped[ni, nj]) base float(np.mean(vals)) # 已展开邻域的平均相位 k round((base - phi[i, j]) / (2.0 * np.pi)) # 最近整数个 2π unwrapped[i, j] phi[i, j] 2.0 * np.pi * k for di, dj in ((1, 0), (-1, 0), (0, 1), (0, -1)): ni, nj i di, j dj if 0 ni hgt and 0 nj wdt and not visited[ni, nj]: heapq.heappush(heap, (-q[ni, nj], ni, nj)) return unwrapped优先队列里用负质量当键pop 出来就是质量最高的待展开像素。每展开一个像素用已展开邻域相位均值做基准计算当前包裹相位与基准相差的 2π 整数倍并补上这就是展开的本质。这个实现里没有处理遮挡和断崖真实工程还要加掩膜与残差点标记但作为仿真实验的基准已经够用。跑完把 unwrapped 画成伪彩图物体区域颜色应该连续且平滑如果出现一条条接缝说明有跳点没补上。4.2 相位差换算高度等效波长模型与真值对比解包裹得到的绝对相位里还叠着载频项要还原高度得先去掉它。仿真里最干净的做法是同样跑一遍零高度参考平面两幅展开相位相减Δφ(x,y) φ_obj(x,y) - φ_ref(x,y)h(x,y) Λ·Δφ(x,y) / (2π)。这个线性关系只在等效波长模型成立时准确。真实系统里投影仪和相机光轴不平行Δφ 与 h 是接近双曲的非线性关系标定时要拟合高阶多项式仿真实验里没这个必要但你要清楚线性模型是近似。下面把整条流水线封装成一个函数物体和参考平面各跑一遍def run_sim(height_map, seed0): 输入高度场输出展开相位和调制度质量图。 rng np.random.default_rng(seed) phi plane_k * x 2.0 * np.pi * height_map / lam_eq fr [0.5 0.5 * np.cos(phi 2.0 * np.pi * k / N) for k in range(N)] fr np.stack(fr) for k in range(N): fr[k] gaussian_filter(fr[k], sigmablur_sigma) # 参数来自 2.3 节 fr[k] rng.normal(0.0, noise_sigma, (h, w)) fr np.clip(np.round(fr * 255.0), 0, 255).astype(np.uint8) phi_w, B phase_shift_reconstruct(fr) return quality_guided_unwrap(phi_w, B), B phi_obj, B_obj run_sim(h_true, seed1) phi_ref, B_ref run_sim(np.zeros_like(h_true), seed2) h_est (phi_obj - phi_ref) * lam_eq / (2.0 * np.pi)两个坑在这里最容易出现。一是全局 2π 常数解包裹只保证相对连续物体和参考平面的展开相位可能差一个 2π 的整数倍导致整幅高度偏移一个 Λ。仿真里如果发现 h_est 整体偏高或偏低先检查这个而不是去调滤波参数。二是 run_sim 用到的 plane_k、x、blur_sigma 都来自第 2 节定义的全局环境封装函数时不要另起一套参数否则两次仿真的空间基准就对不齐。提示先做一步零高度参考平面残差检查展开相位相减后残差均值应接近 0若接近 ±2π就是全局常数漂移。4.3 用 RMSE 和边界掩膜量化重建误差仿真实验的最后一环也是它区别于真机调参的地方和真值 h_true 逐像素比大小。valid mask (B_obj 0.3) # 物体内部且调制度足够 err h_est - h_true rmse np.sqrt(np.mean(err[valid] ** 2)) print(f有效像素 {valid.sum()}, RMSE {rmse:.4f})掩膜必须加原因有两层物体外面的像素没有真值可对而物体边缘因离焦模糊B 骤降误差会大一个数量级不剔除的话 RMSE 会被边缘主导看不出算法本身的水平。想看误差的空间分布就把 err 用有效掩膜画热力图通常边缘一圈亮、内部接近零如果内部出现周期性亮带那就是谐波泄漏或解包裹跳线和随机噪声的形状完全不同。各处理阶段的误差来源整理如下阶段主要误差来源仿真里怎么看条纹生成离焦、量化、加性噪声对比加噪前后 B 的直方图包裹相位谐波泄漏、N 过小用 N12 的结果当近似真值解包裹梯度超过 π/px、低质量区检查展开相位的接缝高度映射全局 2π 常数漂移参考平面区域残差均值5. 相移轮廓术仿真的收尾技巧固定种子、中间结果与 git 管理5.1 固定随机种子并逐级落盘仿真实验改一个参数就要重跑一次如果随机种子每次不一样两轮结果的差异分不清是参数引起的还是噪声引起的。我习惯在 run_sim 里把 seed 作为显式参数传入做对比时固定 seed 只改被测参数同时把包裹相位、展开相位、高度估计逐级用 np.save 落盘排查时不用把整条链路重跑一遍。建议的文件布局是 data/ 放条纹与相位中间结果out/ 放高度图和误差图脚本根目录只留源码。5.2 解压校验、EOCD 报错与 git 管理拿到压缩包先做完整性校验再动手改代码。unzip -t 会逐文件测试 CRC如果报 could not find EOCD说明文件没下完整直接用下载工具重试别浪费时间找所谓修复工具。遇到 .z01 分卷或报 CRC 错误的文件用 7z x 统一处理更省事。unzip -t 相移轮廓术仿真实验.zip unzip 相移轮廓术仿真实验.zip -d psp_sim cd psp_sim git init git add -A git commit -m psp simulation baseline解压后第一件事就是 git init。仿真的参数一改结果就变没有版本记录的话一两周后你根本说不清哪张误差表对应哪套参数。每次调完参提交一次commit message 里写清楚改了 p、N 还是噪声对比实验就能随时回退。最后一个验证技巧把展开相位减去载频项 plane_k·x得到的残差应该是平滑的物体相位如果残差图上有锯齿状条纹说明某个像素的 2π 补跳补错了优先检查 B 最低的那片区域。本文还有配套的精品资源点击获取