全景渲染与融合（投影、缝合、曝光补偿）（图像拼接系列 · 第 7 篇）

你可以把这一篇理解成：前面几篇在解决“几何对不对”，这一篇在解决“最后看起来好不好”。

前几篇我们已经完成了从“相机模型、标定、畸变矫正、特征匹配”
到“匹配图 + MST 初始化 + BA 全局优化”的整条几何链：

• 每一张图像都有了一个较为一致的全局姿态 $R_i$；
• 相机内参 $K$ 已知，柱面/球面等投影模型也已明确；
• 图与图之间的几何关系在公共坐标系下大体自洽。

最后一步，也是用户最直观感受到的一步，就是：

如何利用这些姿态，把所有输入图像重投影到同一个输出面，
并通过缝合、曝光补偿与后处理，生成一张观感自然的全景图？

本篇会回答：

• 输出投影面如何选择？flat / cylindrical / spherical 有何区别？
• 如何根据所有视角范围估计输出全景的尺寸与分辨率？
• 多张图如何在输出面上重投影并融合？
  （简单 alpha、线性渐变、多频段 / 金字塔融合）
• 如何做曝光/增益补偿，减弱接缝处的亮度/色彩跳变？
• 拼接完如何裁剪黑边、做基本后处理？

如果把这一篇压缩成一条非常工程化的主线，其实就是三步：

1. 先决定：输出坐标系长什么样、画布有多大
2. 再决定：输出每个像素该去哪些输入图里采样
3. 最后决定：多张图采到一起时，怎么融合、怎么补曝光、怎么裁边

所以渲染阶段虽然看起来像“后处理”，但本质上是把前面所有几何结果真正落成可看的图像。

---

0 动机：从“姿态”到“图像”的最后一跳

复盘一下，我们目前拥有的信息：

• 每张输入图像 $I_i$ 的像素尺寸 $(w_i,h_i)$；
• 相机内参 $K$（或等效焦距等参数）；
• 每张图像对应的全局姿态 $R_i$（BA 优化结果）；
• 拼接模式：平面（flat）、柱面（cylindrical）、球面（spherical）。

要输出一张全景图 $P$，我们至少要回答三个问题：

1. 输出坐标系是什么？

   • 是平面透视？
   • 还是柱面展开？
   • 或是球面（等距矩形投影）？

2. 输出图像有多大（宽高、分辨率）？

   • 横向覆盖多少视场角？
   • 纵向覆盖多少视场角？
   • 每度多少像素（即采样密度）？

3. 对输出图上每一个像素，该从哪几张输入图的哪个像素位置采样？

   • 如果多个输入图覆盖同一输出像素，如何进行融合？

本篇就是围绕这三点展开的。

---

1 输出投影面与尺寸估计

1.0 三种主流输出投影

与第 3 篇中提到的“投影模式”保持一致，全景输出常用三种：

模式	名称	简介
flat	平面透视	类似单张大广角照片，适合小视场
cylindrical	柱面投影	把场景投影到圆柱再展开，适合水平视场较大、垂直视场中等的全景
spherical	球面投影（equirectangular）	把场景投影到球面再展开，适合 360° 全景或上下覆盖较大的情况

在常见的全景拼接实现中：

• 摄影常见为手持环拍，水平 FOV 很大，垂直 FOV 中等；
• 柱面投影常作为默认输出方式。

---

1.1 以方向向量为中介

无论哪种输出投影，本质上都是：

把每个相机视图中的像素，
通过相机姿态 $R_i$，映射到一个全局方向向量 $\mathbf{v}$ 上，
再把这些方向映射到目标投影面坐标。

在第 6 篇中我们已经用了“方向一致性”的视角，
这里继续沿用这一中间表示。

对图像 $i$ 中的像素 $p_i=(u,v)$：

1. 通过内参和（如果需要的话）投影逆变换，得到相机坐标系下方向：

   $$d_i=\text{PixelToDirection}(p_i,K,\text{输入投影模式})$$

2. 通过姿态 $R_i$ 旋转到全局坐标系：

   $$\mathbf{v}=R_i{d}_i$$

全局方向 $\mathbf{v}$ 再被映射到输出投影面，生成输出像素坐标。

---

1.2 输出尺寸估计（以柱面为例）

以柱面投影为例，坐标关系可写为（略去半径尺度）：

$$\theta =atan2(v_x,v_z),\varphi =\arcsin (v_y)$$

其中：

• $\mathbf{v}=(v_x,v_y,v_z)$ 是全局方向向量
• $\theta$ 是水平角度（绕竖直轴）
• $\varphi$ 可看作垂直方向的“高度角”

若我们只做水平方向的大范围全景，常用的简化是：

• 只用 $\theta$ 控制横坐标
• 用 $\varphi$ 或其某种函数（如 $y$ 在柱面坐标）控制纵坐标

为了估计输出宽高，可以：

1. 对所有图像的所有有效像素（或采样像素）：

   • 按前面的方法计算对应的 $\theta ,\varphi$

2. 统计：

   $${\theta }_{\min },{\theta }_{\max },{\varphi }_{\min },{\varphi }_{\max }$$

3. 设定采样密度，例如：

   • 水平方向每度多少像素：$s_{\theta }$（px / deg）
   • 垂直方向每度多少像素：$s_{\varphi }$

4. 得到输出图宽高：

   width  ≈ (θ_max - θ_min) * s_θ
   height ≈ (φ_max - φ_min) * s_φ
   

工程上通常不会真的遍历所有像素，而是：

• 只在特征点、网格采样点或图像边缘进行采样
• 再通过一定安全裕度（margin）扩展一点范围

一个很重要的经验是：尺寸估计宁可略大、后面再裁，也不要一开始估得过小。

因为：

• 估小了会直接截掉本该保留的画面内容
• 估大了最多只是多出黑边，而黑边后面还可以靠掩膜裁剪掉

所以在工程上，EstimateViewRange 通常追求的是“稳妥覆盖”，而不是数学上刚好卡边。

---

1.3 flat / spherical 情形

• flat（平面透视）：

  • 适合整体视场较小、或只截取其中一段的情况
  • 可以选一个参考视角 $R_{\text{ref}}$ 作为“相机方向”，
    以其为中心定义输出平面；
  • 再根据希望的 FOV 与像素密度计算宽高。

• spherical（球面等距矩形）：

  • 直接使用：

    $$\theta =atan2(v_x,v_z),\varphi =\arcsin (v_y)$$

  • $\theta \in [-\pi ,\pi ]$，$\varphi \in [-\frac{\pi }{2},\frac{\pi }{2}]$

  • 输出宽高可以直接设置为：

    width  = s * 2π
    height = s * π
    

    其中 $s$ 为每弧度多少像素。

---

2 从全局方向到输出像素：几何链路

2.0 输出像素反向推回全局方向

在第 3 篇中我们讲过“反向映射 + 插值”的思路：
在渲染输出图时，常用的是：

对输出每一个像素，反推它对应的全局方向，再查回输入图的采样位置。

原因同之前：

• 输出图每个像素都能确保有且仅有一次赋值
• 不会出现空洞，插值操作也更统一

以柱面输出为例，对输出图的像素坐标 $(x,y)$：

1. 先将其归一化到“柱面参数空间” $(\theta ,\varphi )$：

   θ = (x / width)  * (θ_max - θ_min) + θ_min
   φ = (y / height) * (φ_max - φ_min) + φ_min
   

2. 再把 $(\theta ,\varphi )$ 转换为全局方向 $\mathbf{v}$，
   例如简单球面近似：

   $$\mathbf{v}=\left[\begin{array}{c}\cos \varphi \sin \theta \\ \sin \varphi \\ \cos \varphi \cos \theta \end{array}\right]$$

3. 以 $\mathbf{v}$ 为中间变量，再去查各个输入图的采样位置。

---

2.1 从全局方向到某一图像的像素

对于输入图像 $i$，给定全局方向 $\mathbf{v}$：

1. 转回该相机的局部坐标系：

   $$d_i=R_i^{-1}\mathbf{v}$$

2. 通过针孔模型 / 柱面逆投影等，得到该图像上的像素坐标：

   (u, v) = DirectionToPixel(d_i, K, 输入模式)
   

3. 检查 $(u,v)$ 是否落在图像 $i$ 的有效区域内：

   • 若越界，则图像 $i$ 不覆盖该方向
   • 若在边界附近，可应用一定的外扩/裁剪策略

4. 若覆盖，则可以在 $(u,v)$ 处进行双线性（或更高阶）插值，得到该图的采样颜色。

---

2.2 多图覆盖：可见性与权重

对于输出图某个像素对应的方向 $\mathbf{v}$，
通常会有多张图 $I_i$ 覆盖该方向：

• 比如环拍中，中间一段墙会被多帧重复拍到

此时我们需要：

1. 确定哪些图像 $i$ 有效覆盖该方向：

   • 即 DirectionToPixel 后的 $(u,v)$ 在图像有效区域

2. 为每个有效图像分配一个融合权重 $w_i$：

   • 可以与视角距离、视线角度、接缝设计等因素有关

3. 按照权重对采样颜色做加权平均：

   $$C_{\text{out}}=\frac{{\sum }_i{w}_i{C}_i}{{\sum }_i{w}_i}$$

这里的 $C_i$ 是图像 $i$ 在 $(u_i,v_i)$ 位置的颜色（插值结果）。

一个很常见的工程写法，是把权重拆成几个可解释的因子：

$$w_i=w_i^{\text{border}}\cdot w_i^{\text{angle}}\cdot w_i^{\text{seam}}$$

例如：

• $w_i^{\text{border}}$：越靠近图像中心越大，越靠近边缘越小
• $w_i^{\text{angle}}$：该方向越接近相机主视轴，权重越大
• $w_i^{\text{seam}}$：若某处被 seam 选择为“尽量不要从这张图取”，则权重变小

这样做的好处是：每个因子都对应一个明确的工程意图，后续调参也更清楚。

把这一节压成一句话，就是：

渲染的几何主线是：输出像素 → 全局方向 → 各输入图像像素 → 采样颜色 → 按权重融合。

---

3 融合策略：从简单 alpha 到多频段

3.0 简单覆盖与 alpha 融合

最原始的策略是：

按某个顺序依次把输入图“画”到输出图上，后画的覆盖先画的。

这种做法在重叠区域会出现明显接缝和“硬边界”，
但实现极其简单，适合作为 debug 版本。

略好一点的方式是 alpha 融合：

• 对每个输入图定义一张权重图 $W_i(x,y)$，
  通常在重叠区域中间权重大、边缘权重小

• 最终颜色为：

  $$C_{\text{out}}=\frac{{\sum }_i{W}_i{C}_i}{{\sum }_i{W}_i}$$

若 $W_i$ 在各自图像的左右边缘做线性渐变，
就能在接缝处形成一个过渡带，减小“硬接缝”。

这里顺手区分两个很容易混淆的概念：

• seam finding（找缝）：决定“重叠区域里到底从哪张图为主”
• blending（融合）：决定“在选定的缝附近，怎样把两张图过渡得更自然”

很多时候画面难看，不是因为“不会融合”，而是因为缝选错了地方——比如正好切过一条强边缘、一扇窗框或一个人脸。

所以比较成熟的流水线往往是：

1. 先用代价函数/图割等方法找一条相对不显眼的 seam
2. 再在 seam 附近做 feathering 或 multi-band blending

---

3.1 线性渐变缝合（feathering）

典型 feathering 策略：

• 对每张投影后的输入图 $I_i$，
  为其定义一个“距离图像边界的距离” $d_i(x,y)$

• 权重取为：

  $$W_i(x,y)=\min \left(1,\frac{d_i(x,y)}{D}\right)$$

  其中 $D$ 为预设的衰减宽度

直觉：

• 在图像中心 $d_i$ 大，权重趋近 1
• 在边缘的 $D$ 像素范围内，权重从 0 线性过渡到 1
• 重叠区域中，两张图的过渡带会叠加，形成平滑缝合

伪代码示例（已经假设输入图都投影到同一输出坐标系）：

for each input image i:
    for each pixel (x, y) in its bounding box on panorama:
        if pixel valid:
            d = DistanceToImageBoundary(x, y, image_i)
            w = clamp(d / D, 0, 1)
            C = SampleColor(image_i, x, y)

            acc_color[x, y] += w * C
            acc_weight[x, y] += w

for each output pixel (x, y):
    if acc_weight[x, y] > 0:
        pano[x, y] = acc_color[x, y] / acc_weight[x, y]
    else:
        pano[x, y] = background_color

---

3.2 多频段 / 拉普拉斯金字塔融合（multi-band blending）

在光照差较大、场景结构复杂时，简单的 feathering 仍然会露馅：

• 低频部分（大范围亮度差）容易形成可见的宽接缝
• 高频细节（边缘、纹理）在接缝附近可能出现重影或模糊

多频段融合的核心思想是：

把图像分解成多个尺度（频段），
在低频层使用较宽的过渡带、
在高频层使用较窄（甚至锐利）的过渡带，
最后再把各频段融合结果重建回来。

流程大致为：

1. 对每张输入图构建高斯金字塔 $G_i^0,\dots ,G_i^L$
2. 由高斯金字塔构建拉普拉斯金字塔 $L_i^0,\dots ,L_i^L$
3. 对每一层 $l$，按某种权重图 $W_i^l$ 对 $L_i^l$ 做加权融合
4. 从最高层开始逐级向下重建，得到最终融合图

本文不展开完整公式，只保留直觉：
“低频负责整体亮度/色彩的一致过渡，高频负责边缘/细节的锐利融合”。

在实际拼接库中，多频段融合往往是提升观感最明显的一步，
尤其是在室外场景（天空 + 建筑）、室内灯光复杂时。

---

4 曝光与增益补偿

即使几何对齐与融合策略都做得很好，
你仍然会看到这样的问题：

• 不同照片之间的曝光 / ISO / 白平衡略有变化
• 重叠区域整体偏亮 / 偏暗 / 偏黄 / 偏蓝
• 结果是：接缝处出现明显的“颜色带”

曝光补偿 / 增益估计的目标就是：

在保证整体动态范围的前提下，
估计每张输入图的一个简单亮度/颜色变换（如增益、偏置），
使得重叠区域的统计量尽量一致。

---

4.0 灰度增益的最简模型

从最简单的灰度增益模型开始：

$$I_i^{\prime }(x,y)=g_i\cdot I_i(x,y)$$

其中 $g_i$ 为图像 $i$ 的亮度增益（实数）。

在两张图 $i,j$ 的重叠区域中，我们希望：

$$g_i{I}_i\approx g_j{I}_j$$

可以用最小二乘意义下的约束来估计所有 $g_i$。

---

4.1 在线性方程组中估计增益

对重叠区域中采样的若干像素 $(x_k,y_k)$，
我们有：

$$g_i{I}_i(x_k,y_k)\approx g_j{I}_j(x_k,y_k)$$

可改写为：

$$g_i{I}_i(x_k,y_k)-g_j{I}_j(x_k,y_k)\approx 0$$

对所有重叠像素、所有相邻图像对 $(i,j)$，
可以累加出一个线性最小二乘问题，
未知量是所有 $g_i$。

由于整体放大/缩小对结果影响不大，
通常会增加一个约束：

$$\sum _i{g}_i=N\text{或者}{g}_{\text{ref}}=1$$

其中 $N$ 为图像数量，ref 为参考图。

求解得到 $g_i$ 后，在渲染/融合前对每张图做一次增益校正即可。

如果你更希望把问题写成“标准线性”的形式，也常见一种对数域写法。令：

$$a_i=\log g_i$$

则在像素强度为正时，约束

$$g_i{I}_i\approx g_j{I}_j$$

可改写为近似的：

$$a_i-a_j\approx \log I_j-\log I_i$$

这样未知量就从乘法增益变成了加法量，更便于直接塞进线性最小二乘框架。

---

4.2 RGB 通道与颜色偏差

若要处理颜色偏差，可以扩展为每通道一个增益：

Ii′(x,y,c)=gi,c⋅Ii(x,y,c),c∈{R,G,B} I_i'(x,y,c) = g_{i,c} \cdot I_i(x,y,c), \quad c \in \{R,G,B\} Ii′​(x,y,c)=gi,c​⋅Ii​(x,y,c),c∈{R,G,B}

构造的线性系统与灰度情形类似，只是变量数量变为 $3N$；
也可以在 YUV/HSV 等颜色空间中仅对亮度通道做增益，
对色度通道做较弱的调整。

实践中常见折中：

• 先在亮度通道上做主增益补偿（修掉亮度带）
• 再对色彩做较小幅度的平滑/白平衡调整

---

4.3 全局 VS 局部增益

上述方法是per-image 全局增益，
在某些极端情况下（例如强烈的曝光变化、局部阴影）仍然不够。

可选的更复杂方案是：

• 对图像划分网格（例如 $m\times n$ 区块）
• 为每个区块估计一个增益 $g_{i,block}$
• 再在图像内部对这些增益做平滑插值

这类方法虽然更灵活，但也更易过拟合和产生伪影，
需要更复杂的正则化设计，这里不展开。

---

5 边界裁剪与后处理

5.0 黑边与空洞区域

无论是柱面还是球面投影，
在重投影之后，你常会看到：

• 上下或四周有大片黑边（无像素覆盖区域）
• 部分区域被极度拉伸，信息贫乏

典型处理方式：

1. 有效区域掩膜：
   在渲染时同步维护一张“有效像素掩膜”，
   标记哪些输出像素至少被一张输入图覆盖。

2. 边界裁剪：

   • 对掩膜取连通区域，找到最大连通块
   • 对该连通块求其外接矩形 / 凸包
   • 根据用户需求裁剪到适当矩形范围

这样可以自动去掉大部分纯黑 / 伪影区域。

---

5.1 简单的色调映射与锐化

完成几何、融合与曝光补偿后，
还可以做一些轻量的后处理：

• 全局色调映射（tone mapping）：
  • 拉伸对比度
  • 适度增加饱和度
• 轻微锐化：
  • 如 Unsharp Mask 等
  • 提升细节观感

注意不要过度：

• 过强的锐化会放大噪声与接缝残留
• 过激的色调映射会产生 halo 与色彩断层

---

6 典型渲染流水线伪代码

为了把前几节串起来，这里给出一个简化版的整体伪代码。

function RenderPanorama(images, K, R_final, config):
    # 1. 估计输出投影范围与尺寸
    θ_min, θ_max, φ_min, φ_max = EstimateViewRange(images, K, R_final, config.projection_type)
    (W, H) = ComputePanoramaSize(θ_min, θ_max, φ_min, φ_max, config.resolution)

    pano_color   = zeros(H, W, 3)
    pano_weight  = zeros(H, W)
    pano_valid   = zeros(H, W, bool)

    # 2. 预先估计每张图的曝光增益 g_i
    gains = EstimateExposureGains(images, K, R_final, config)

    # 3. 主循环：对输出每个像素反向查找
    for y in 0..H-1:
        for x in 0..W-1:
            (θ, φ) = PixelToGlobalAngles(x, y, θ_min, θ_max, φ_min, φ_max, W, H, config.projection_type)
            v = AnglesToDirection(θ, φ, config.projection_type)  # 全局方向

            blended_color = 0
            blended_weight = 0

            for each image i:
                R_i = R_final[i]
                d_i = inverse(R_i) * v

                (u, v_img) = DirectionToPixel(d_i, K, config.input_projection_type)
                if not InsideImage(u, v_img, images[i].size):
                    continue

                C = BilinearSample(images[i], u, v_img)
                C = gains[i] * C    # 曝光增益

                w_i = ComputeBlendWeight(i, x, y, u, v_img, config)

                blended_color  += w_i * C
                blended_weight += w_i

            if blended_weight > 0:
                pano_color[y, x] = blended_color / blended_weight
                pano_weight[y, x] = blended_weight
                pano_valid[y, x]  = true

    # 4. 裁剪黑边 / 无效区域
    crop_rect = ComputeCropRectFromMask(pano_valid, config)
    pano_color = Crop(pano_color, crop_rect)

    # 5. 可选后处理（tone mapping, sharpening）
    pano_color = PostProcess(pano_color, config)

    return pano_color

实际实现中，多数操作会搬到 GPU 上，
并使用更复杂的融合与增益估计策略，上述伪代码强调的是几何与数据流主线。

---

7 常见坑与 FAQ（渲染与融合阶段）

---

7.1 症状：接缝处几何对齐很好，但亮度 / 颜色明显不一致

• 常见原因：
  • 曝光/白平衡差异没有做任何补偿
  • 只做了简单的 alpha 混合，未在重叠区域调整统计量
• 建议：
  • 至少对重叠区域做 per-image 灰度增益估计
  • 在高差异区域使用更宽的融合带或多频段融合

---

7.2 症状：接缝处几何上有重影或模糊

• 常见原因：
  • BA 后某些局部姿态仍有残差（几何上没完全对齐）
  • 融合带太宽，导致两边略错位的结构被平均
  • 使用简单线性混合处理复杂几何（如近景、明显视差）
• 建议：
  • 检查 BA 残差分布，排查局部不收敛的视图
  • 对结构复杂区域减少融合带宽度，甚至使用“硬切 + 细致 seam carving”方案
  • 对近景/视差极大的数据，本身就超出“纯旋转拼接”的适用范围

---

7.3 症状：输出图上下（或左右）存在大块空白或严重拉伸

• 常见原因：
  • 输出角度范围估计过宽，包含了几乎无数据覆盖的方向
  • 采样密度设置不合理（某方向过密，导致拉伸）
• 建议：
  • 在 EstimateViewRange 阶段更严格地仅考虑有效方向
  • 渲染后通过掩膜自动裁剪；
  • 或在 UI 中允许用户交互式裁剪

---

7.4 症状：整体色调怪异或对比度异常

• 常见原因：
  • 曝光增益估计过度，导致整体亮度被扭曲
  • 后处理（tone mapping / sharpening）参数过激
• 建议：
  • 在引入复杂增益/颜色模型之前，先用最简单的灰度增益调通主链路
  • 在小样本上可视化增益分布，避免极端值
  • 逐步调整后处理强度，观察变化趋势而非一次性调大

---

8 本系列小结与展望

这一篇是本系列的“收尾篇”，我们终于从理论和工程角度
走完了一条典型的全景拼接管线：

1. 相机建模与标定
   • 针孔模型、内参 $K$、外参 $(R,t)$，
   • 畸变模型 $D$ 与棋盘标定。
2. 畸变矫正与投影
   • 归一化坐标、去畸变、反向映射 + 插值、
   • 柱面/球面投影与 half-shifted 坐标。
3. 特征匹配与单应估计
   • 特征点与描述子、ratio test、
   • DLT 估计单应 $H$ 与 RANSAC 稳健拟合、
   • 在纯旋转假设下从 $H$ 求 $R$。
4. 匹配图与 MST 初始化
   • 构建 view graph，
   • 以边质量为代价构造 MST，
   • 沿树传播相对旋转 $R_{ij}$ 得到初始 $R_i$。
5. BA 全局优化
   • 把所有匹配约束变成“方向一致性”的残差，
   • 在旋转向量参数化下用 Gauss-Newton / LM 做全局优化。
6. 渲染与融合（本篇）
   • 选择输出投影面与分辨率，
   • 方向 → 像素的正反几何链路，
   • 多图重投影、融合与曝光补偿，
   • 边界裁剪与后处理。

如果你能跟下来并把伪代码翻译成自己的实现，
基本就能从零搭建一个“可用级别”的全景拼接系统。

---

8.1 自测（全系列回顾）

1. 单应矩阵与纯旋转
   在什么条件下，两张图像之间的关系可以用 $H=KR{K}^{-1}$ 近似？你会如何从 $H$ 反推出 $R$？
2. 反向映射与 LUT
   为什么工程中更偏爱“反向映射 + 插值”？LUT 在其中起什么作用？
3. RANSAC 的迭代次数公式
   回忆内点率 $w$、采样大小 $s$ 时，成功概率 $P$ 与迭代次数 $N$ 的关系。
4. MST 初始化 vs BA 优化
   用自己的话说明：MST 初始化给了什么，BA 又在此基础上做了什么？
5. 渲染融合中的几何链
   从输出像素 $(x,y)$ 出发，写出“输出像素 → 全局方向 → 各图像像素 → 采样与融合”的完整几何链条。

---

8.2 符号速查（本篇新增）

符号	含义
$R_i$	图像 $i$ 的全局旋转姿态（BA 后结果）
$\mathbf{v}$	全局方向向量（单位向量）
$(\theta ,\varphi )$	方向 $\mathbf{v}$ 在柱面/球面投影下的角度坐标
$W,H$	输出全景图的宽和高
$W_i(x,y)$	输入图 $i$ 在输出坐标 $(x,y)$ 处的融合权重
$g_i$	输入图 $i$ 的曝光增益
$C_i$	输入图 $i$ 在某个采样位置的颜色向量
$C_{\text{out}}$	输出全景在某像素位置的最终颜色