简介：这是一套基于光度立体技术的三维重建应用程序Python源码包，面向计算机视觉方向的课程设计、毕业设计及入门进阶学习者。项目从多角度光照图像出发，完成法向量估计、深度图恢复与三维模型生成，包含完整可运行代码与详细注释，配套项目说明、光度立体算法流程图、实验数据和结果展示文件，适合直接运行验证，也支持二次开发。压缩包共41个文件，以png结果图、py脚本、ipynb交互式笔记本、npy数据文件、obj三维模型及txt说明文档为主，另含csv/xls辅助数据，整体仅7.67MB，结构清晰、便于按流程学习。目前已有243人学习下载。使用时可先阅读“项目必看.txt”，将解压目录改为英文路径，再按main.ipynb逐步执行，即可复现从图像输入到三维重建输出的完整流程。

1. 光度立体三维重建：12张猫图恢复法向量与深度的Python源码

光度立体（Photometric Stereo）是从不同光照方向的多张图像中恢复物体表面法向量与深度的技术，核心诱惑在于：不需要双目相机、不需要结构光、不需要深度传感器，一台固定相机加多个光源就能逐像素算出表面几何。这套Python源码正是围绕这个思路做的完整实现，输入是同一只猫在12个不同光照方向下的图片，输出是法向量图、深度图和cat.obj模型。你要是正在做计算机视觉课程设计、毕业设计，或者想亲手把光度立体算法跑通再二次开发，这份源码可以直接复现。整个工程带详细中文注释，数据、脚本、notebook、说明文档齐全，还附带一张算法流程图给答辩PPT用。

2. 光度立体的原理与项目文件结构：光照矩阵和掩码决定重建上限

在跑代码之前，先花十分钟弄明白算法在算什么。光度立体不是SLAM那种多视角几何路线，它只依赖一个物理模型：被拍摄物体表面近似为朗伯体，也就是漫反射表面，表面亮度只与入射光方向和表面法向量的夹角有关。理解了这一点，后面调参数时你就知道哪些问题该查数据，哪些问题该查代码。

2.1 朗伯体假设与光度学方程：为什么12张图比3张图稳

对于固定机位下的同一个物体，第i张图像在像素(x,y)处的灰度可以写成：

I_i(x,y) = ρ(x,y) · N(x,y) · L_i

其中ρ是表面反照率（albedo），N是待求的表面单位法向量，L_i是第i张图像对应的光照方向单位向量。把12个光照方向堆成一个3×12的光照矩阵，再对每个前景像素收集12个观测灰度值，就得到一个超定线性方程组。按最小二乘或伪逆求解，先解出带反照率的ρN，再取模长分离出ρ和单位法向量N。

从数学上看，三张非共面的光照方向就能让矩阵满秩，但实际拍摄中单张图像必然有传感器噪声、暗电流和阴影边缘误差。我一般建议至少拍8到12张，光源方向尽量布在半球空间，不要把所有光源挤在一个小范围内，否则光照矩阵的条件数会变大，最小二乘解对噪声更敏感。这套猫数据正好是12张，冗余足够，又不至于让拍摄和计算负担过重。还有一个隐含前提：整个算法建立在正交投影或弱透视假设上，相机视角不能太广，物体深度变化不能太大。所以拍摄时相机要拉远一点，让物体占画面大约二分之一到三分之二，贴近物体拍的广角图会让重建结果出现系统性畸变。

2.2 项目资源结构：每个文件在流程里扮演什么角色

压缩包解压后，核心模块是Photometric_Stereo.py，它是光度立体重建的主程序，负责读取图片、解析光照矩阵、求解法向量、积分深度并导出OBJ。iamge_operations.py是图像预处理辅助模块，提供灰度化、掩码处理、归一化这类基础操作。main.ipynb是逐步骤可交互运行的notebook，适合边看结果边理解算法，也方便做课程演示。数据都在untreated_data目录下，根目录还有3.jpg、4.jpg这类效果展示图，以及normal_vactor_cat.jpg、depth_cat.jpg、cat.obj这些运行产物。

输入文件与作用的对应关系如下：

目录里还有一张光度立体算法流程图.png，是给论文或答辩PPT准备的结构图素材。压缩包里附带了一个遗传算法优化神经网络的完整项目（GA-NN.ipynb加实验报告），跟三维重建无关，属于额外材料，用到的时候再单独解压。另外你会看到__pycache__下有Photometric_Stereo.cpython-38.pyc，说明这套代码在Python 3.8环境下已经实际运行过，环境兼容性不需要担心。

这里有一个值得养成的习惯：npy文件是重要中间产物，二次开发时不用重新跑一遍12张图的光度立体求解，直接np.load('normal_vactor_cat.npy')就能拿到法向量场。把中间数据单独落盘，调试深度积分算法时能省掉大量重复计算时间，我就是这么干的。

2.3 光照矩阵light.txt怎么读：列顺序决定重建方向是否正确

light.txt是整个流程里最容易出错的文件。它保存的是一个3×N的矩阵，N等于图像数，每一列是某张图像对应的光照方向单位向量，形式是(x, y, z)，其中z轴指向相机。打开文件看到的格式类似下面这样：

0.000000 0.000000 1.000000 0.183013 0.683013 0.707107 -0.500000 0.500000 0.707107

第一列对应cat.0.png，第二列对应cat.1.png，顺序必须和cat.txt里列出的文件名一一对应。判断光照方向是否正确有个经验办法：看重建出的法向量图。把法向量的x、y、z分量分别映射成R、G、B通道，如果猫的鼻子区域显示为偏蓝色，说明法向量朝向z轴也就是朝向相机，额头到鼻尖的颜色过渡平滑，说明光照矩阵方向和顺序都对。如果法向量图出现大面积异常渐变或颜色整体偏移，优先怀疑light.txt的列顺序跟图像名没对上，而不是先去怀疑积分算法。

光源方向是所有像素共享的全局量，它对整张图都生效。如果你的拍摄过程中光源位置在图像之间发生了移动，但光照矩阵没同步更新，重建出的法向量场就会朝某一方向系统性偏移，深度图也会跟着变形。这也是为什么强烈建议把图像名与光源方向的映射关系打包保存，而不是单独丢一个light.txt在目录角落里。

3. 核心代码拆解：Photometric_Stereo.py里的最小二乘与深度积分

这一章把主程序的三个关键阶段逐行拆开。代码结构和我在实际项目里常用的实现基本一致，只是这个工程把图像读取、掩码筛选、法向量估计、深度积分、OBJ导出拆成了几个独立函数，方便单步调试。以下代码是结合项目主流程梳理出的核心写法，实际文件里注释更密。

3.1 图像读取与掩码筛选：为什么用float64而不用uint8

主程序首先会扫描untreated_data目录，按cat.txt里的文件名列表读取所有图片，同时读取mask。iamge_operations.py把这一步封装成了加载函数。

import numpy as np from PIL import Image def load_source_images(img_paths, mask_path): # img_paths: 按cat.txt顺序整理好的图像路径列表 imgs = [] for p in img_paths: # convert('L')转单通道灰度，float64防止运算截断 im = np.asarray(Image.open(p).convert('L'), dtype=np.float64) imgs.append(im) I = np.stack(imgs, axis=-1) # 形状 (H, W, N) mask = np.asarray(Image.open(mask_path).convert('L'), dtype=np.float64) mask = (mask > 128).astype(np.float64) return I, mask

这里有两个容易被忽略的参数选择。第一个是convert('L')把图像转成单通道灰度，光度立体只需要亮度信息，彩色图不参与计算。第二个是dtype用float64而不是uint8，因为后续要做矩阵除法和减法，uint8在计算过程中一旦出现超过255或小于0的值就会截断，产生大量环状伪影。mask的处理方式是直接以128为阈值二值化，确保mask只含0和1两个值。

调用时常见做法是先把cat.txt读成列表，再拼上目录前缀：

import os img_dir = 'untreated_data' with open(os.path.join(img_dir, 'cat.txt'), 'r') as f: img_names = [line.strip() for line in f if line.strip()] img_paths = [os.path.join(img_dir, n) for n in img_names] I, mask = load_source_images(img_paths, os.path.join(img_dir, 'cat.mask.png')) print(I.shape) # 期望输出 (H, W, 12)

如果print出的第三维不是12，说明cat.txt列出的图像数和light.txt的列数不一致，这会在后面的lstsq阶段直接抛维度错误。所以先把三者的维度关系对齐，是运行这套代码的第一道检查线。

3.2 最小二乘求解法向量：lstsq的维度变化与反照率分离

核心求解段在每个前景像素上做一次最小二乘。对像素(x,y)，12张图的灰度值组成观测向量b（长度12），光照矩阵的转置构成系数矩阵A（形状12×3），要解的是三元向量m = ρN。

def estimate_normals(I, mask, L): H, W, N = I.shape L = np.asarray(L).T # (N, 3) normals = np.zeros((H, W, 3)) albedo = np.zeros((H, W)) ys, xs = np.where(mask > 0.5) for x, y in zip(xs, ys): b = I[y, x] # (N,) m, _, _, _ = np.linalg.lstsq(L, b, rcond=None) rho = np.linalg.norm(m) # 反照率，即m的模长 if rho > 1e-8: normals[y, x] = m / rho albedo[y, x] = rho return normals, albedo

说明一下参数：rcond=None是NumPy 1.14以后的推荐写法，它让lstsq用机器精度去截断奇异值，防止光照矩阵接近奇异时解出巨大数值。rho是m的二范数，物理含义就是表面反照率。除以rho之后把法向量归一化到单位长度。rho > 1e-8这个阈值是防止黑色像素除零的兜底，小于阈值的像素直接保持法向量为零向量，这类像素通常是mask内部反照率极低的黑点。

这版循环实现的好处是每一行每一列都可以单独打印、单独调试，坏处是像素多了之后速度起不来。如果图像是500×500并且mask覆盖一半像素，就是十几万次lstsq调用，跑下来要几秒到几十秒。实际项目中我会把它向量化，这里提一下思路：把mask内的所有像素灰度先取出来，拼成一个大矩阵b_all（形状N×P），再一次性np.linalg.pinv(L) @ b_all得到所有像素的解，速度比循环快一个量级。首次跑通建议用循环版，确认结果正确后再换向量化版本提速。

3.3 深度恢复与OBJ导出：从法向量场积分出高度场

拿到法向量场之后，还需要把方向信息变成z(x,y)高度。这里用的是梯度积分思路：表面梯度可以由法向量直接换算，再从梯度重建高度。

def normal_to_depth(normals, mask): H, W, _ = normals.shape nz = normals[..., 2] + 1e-9 dzdx = -normals[..., 0] / nz # x方向梯度 dzdy = -normals[..., 1] / nz # y方向梯度 z = np.zeros((H, W)) for y in range(H): for x in range(W): if mask[y, x] == 0: continue if x > 0 and mask[y, x-1]: z[y, x] = z[y, x-1] + dzdx[y, x] elif y > 0 and mask[y-1, x]: z[y, x] = z[y-1, x] + dzdy[y, x] z -= z[mask > 0].mean() return z

这个循环积分是最朴素的行扫描方案，只考虑左边和上边已知的深度，速度尚可，但误差会沿扫描方向累积。项目中输出的depth图仔细看的话，可能会发现某一方向上有轻微漂移，这就是行扫描积分的固有特性。更稳的做法是Frankot-Chellappa频域积分，把积分问题转换到傅里叶域，一次性最小化全场梯度误差。我在实际项目里通常先跑一次行扫描积分看整体形态，再用频域积分得到平滑版本。

导出OBJ的代码逻辑同样清晰：遍历mask有效像素，记录顶点坐标，再按四邻域构成三角面。

def export_obj(z, mask, obj_path): H, W = z.shape index_map = np.full((H, W), -1, dtype=np.int32) verts = [] for y in range(H): for x in range(W): if mask[y, x] > 0: index_map[y, x] = len(verts) verts.append((x, y, z[y, x])) with open(obj_path, 'w') as f: for v in verts: f.write("v %.6f %.6f %.6f\n" % v) for y in range(H - 1): for x in range(W - 1): i0 = index_map[y, x] i1 = index_map[y, x + 1] i2 = index_map[y + 1, x] i3 = index_map[y + 1, x + 1] if i0 >= 0 and i1 >= 0 and i2 >= 0: # OBJ的f行索引从1开始，要加1 f.write("f %d %d %d\n" % (i0 + 1, i1 + 1, i2 + 1)) if i1 >= 0 and i3 >= 0 and i2 >= 0: f.write("f %d %d %d\n" % (i1 + 1, i3 + 1, i2 + 1))

index_map记录的是像素坐标到顶点编号的映射，OBJ格式里f行的索引从1开始，所以写入时全部加1。每个四像素邻域拆成两个三角形，共用一个对角线，形成完整的三角网格面。导出后在MeshLab里打开cat.obj，你看到的那个立体猫，就是这12张照片重建出来的几何结果。

4. 常见问题与避坑：路径中文、光照顺序和噪声的排查记录

这一章写的都是我实际跑光度立体项目时踩过的坑，按现象、原因、解决三步整理。每一条都会在高频重复出现，对照着查能省不少时间。

4.1 现象：项目路径带中文，运行报编码错误或FileNotFoundError

这种报错看起来有点玄学，明明文件名是对的，但程序就是找不到文件。原因是Windows控制台的编码和Python源码里的字符串编码不一致，中文字符路径在PIL、numpy读取时偶尔能过、偶尔失败，表现非常随机，跟具体的系统和Python版本都有关系。

解决方式其实很简单：解压后把整个项目改成英文路径，比如D:/photometric-stereo-cat，同时确认文件夹名字也不含中文。项目说明里特别强调过这一点，属于经典的环境坑，不是算法问题。我见过不少人卡在第一步，代码逻辑完全没问题，就是路径里的中文在作祟。

4.2 现象：重建出的深度图出现横向条纹或水波纹

现象是深度图上有一道一道的条纹，像水面波纹一样，而不是平滑的曲面。原因通常是行扫描积分逐行累积误差，任何一行出现一个异常梯度，后续像素全被污染。另一个常见原因是光源方向标定不够准，导致梯度场里混入了低频误差分量。

解决方法是先确认光照矩阵没问题，再把积分算法换成Frankot-Chellappa频域积分。如果不想动算法，可以先对法向量场做一次3×3中值滤波，把明显的异常梯度噪点压下去，再看条纹是否缓解。这个处理对边缘清晰、表面光滑的物体尤其有效。

4.3 现象：法向量图整体偏色，猫鼻子方向明显不对

这是翻车概率最高的环节。现象是法向量图的颜色分布跟预期的R、G、B映射对不上，鼻子朝向不对，整体色偏严重。原因基本可以锁定是light.txt的列顺序和图像顺序错位。光照方向和图像是一一映射的，错一位就相当于在一个像素上叠了错误的光照，最小二乘解出来的不再是真实法向量。

解决方法是写一个自检脚本，逐列打印文件名和对应的光照向量，人工过一遍。我通常还会做一个更直观的检查：把第一行光照向量设为(0,0,1)也就是正对相机方向，对应的那张图应该是最亮的正面光。如果发现cat.0.png不是最亮的那张，那光照顺序肯定有问题。

4.4 现象：mask边缘出现一圈黑白噪点

现象很典型：重建结果的边缘轮廓上，有一圈密密麻麻的黑白像素点，像毛刺一样。原因是mask是硬二值边界，但实际照片里物体边缘存在半影区，像素灰度介于背景和表面之间。这些像素参与最小二乘时，会被解出模长很大且方向混乱的法向量。

解决方法是把mask向内腐蚀几个像素，或者对反照率加置信度权重。项目中可以用scipy.ndimage.binary_erosion实现，腐蚀的半径根据图像分辨率来定，一般1到3个像素就够。腐蚀掉半影区之后，边缘噪点会明显减少，代价是会损失一点真实的边缘几何。

4.5 现象：黑色物体区域重建后法向量大面积为零

比如猫的眼睛、深色花纹，这些区域重建后是一块空洞。原因是黑色表面的灰度接近0，观测向量里的信噪比非常低，lstsq解出的m模长接近0，被rho > 1e-8的阈值过滤掉了。这是朗伯体假设的物理边界：黑色物体几乎不反射光，相机接收不到有效信号。

解决方式有两种。一是提高光源强度、延长曝光时间，让暗部灰度从个位数提升到几十以上；二是干脆在mask里把这类区域去掉，只在能信的区域做重建。我做课程项目时会选择第二种，不做无谓的和物理规律对抗的事情。

5. 进阶：把猫换成自己的物体时的光源标定与结果验证

算法跑通之后，下一步基本就是换自己的数据。这里讲一个我自己验证过的处理顺序。

5.1 自拍一组光度立体数据的现场布置

手机固定在三脚架上，物体放在转台上，保持相机与物体在整个拍摄过程中完全不动。用手电筒或台灯绕物体在半球面上取12个方向，每次只移动光源，不要碰手机和物体。拍完后按顺序编号，手工用画图工具抠出一个二值mask，把背景涂黑、物体涂白。这里有个细节：房间里的环境光要尽量关掉，不然每张图的背景亮度不一致，会影响最小二乘解。

5.2 光源方向标定：没有light.txt时的高效做法

先放一个镜面小球在同一位置拍一张参考图，高光点的位置就是光源方向的指向。把高光点的像素坐标转成方向向量，公式其实是球面几何投影，用球的半径和圆心坐标做反算。不想标定也可以做一个简化验证：规则摆放光源阵时，直接用几何关系算方向，然后观察法向量图的颜色是否平滑。这一步的目的是让light.txt至少有参考价值，而不是纯靠猜。

5.3 用三个指标判断重建结果可不可用

第一个指标是法向量图的平滑度，颜色过渡自然、没有大面积突变，说明光照矩阵基本正确。第二个指标是深度图的灰度直方图，连续分布而不是双峰断裂。第三个指标是把OBJ模型拖进Blender，在材质面板上查看法线朝向，和物体实际的凹凸形状对照。任何一项异常，都优先回到光照顺序和mask边缘检查，而不是改积分算法。

从那以后我每次处理光度立体数据，都强制走一遍先核对光照顺序、再检查mask边缘、最后才看深度结果的流程。这套顺序帮我避开过不少返工，希望这些记录能帮到你。

本文还有配套的精品资源，点击获取