泊松图像融合算法 C++ 实现

泊松图像融合是 Pérez 等人在 2003 年发表的论文 [Poisson Image Editing](
https://www.cs.jhu.edu/~misha/Fall07/Papers/Perez03.pdf)中的算法，该算法实现了两幅图片感兴趣区域的无缝融合。应该是 2015 年 5 月底我在网络上看到这个算法的介绍和结果示意图，当时就对这个算法产生了很大的兴趣。现在网上这个算法的实现太多了，MATLAB C++ java 等语言的实现都有，还收录到了 OpenCV extra module 中。不过我还是自己实现了一遍。

原理

我们可以在目的图片上确定一个感兴趣区域，把这个区域用一张源图片中的区域替代。由于两张图片的差异，直接把源图片中的区域粘贴到目的图片，会出现明显的差异。为此，我们增加一些修改目的图片的约束：第一，目的图片中感兴趣区域的边界像素和外部周边的像素没有显著差异，第二，目的图片中感兴趣区域内部的像素的梯度变化趋势符合源图片的对应像素的梯度变化趋势。见上述论文公式 (7)。这样就达到了无缝内容复制，即 seamless cloning 的目的。

本质上说，该算法是用一组特定的梯度变化图作为引导，在确保边界不变的情况下，修改一张图片某个封闭区域的像素值。

实现

泊松图像融合算法需要求解一个线性方程组 $\mathbf A \mathbf x = \mathbf b$。根据论文的建议，我用 Gauss-Seidel 算法 求解了这个方程。我最早期的实现中，我直接采用 OpenCV 的 cv::Mat 去表示 $\mathbf A$，由于这个矩阵的行列数等于感兴趣区域像素的数量，当感兴趣区域比较大的时候，存储数据量很大，而且运算速度很慢。根据论文中的公式 (7) 很容易发现其实这个矩阵 $\mathbf A$ 是一个稀疏矩阵，于是我自己写了一个 SparseMat 结构，并用这个稀疏矩阵求解方程，代码本文末尾。代码的编译依赖 OpenCV 2.4.5 或者以上版本。

结果示例

示例1

源图片

感兴趣区域

目的图片

融合结果

示例2

源图片

感兴趣区域

目的图片

融合结果

示例3

源图片

目的图片

融合结果

代码

#include <opencv2/core/core.hpp>
#include <opencv2/imgproc/imgproc.hpp>
#include <opencv2/highgui/highgui.hpp>
#include <cmath>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <vector>

struct IndexedValue
{
    IndexedValue() : index(-1), value(0) {}
    IndexedValue(int index_, double value_) : index(index_), value(value_) {}
    int index;
    double value;
};

struct SparseMat
{
    SparseMat() : rows(0), maxCols(0) {}

    SparseMat(int rows_, int cols_) :
        rows(0), maxCols(0)
    {
        create(rows_, cols_);
    }

    void create(int rows_, int cols_)
    {
        CV_Assert(rows_ > 0 && cols_ > 0);

        rows = rows_;
        maxCols = cols_;
        buf.resize(rows * maxCols);
        data = &buf[0];
        memset(data, -1, rows * maxCols * sizeof(IndexedValue));
        count.resize(rows);
        memset(&count[0], 0, rows * sizeof(int));
    }

    void release()
    {
        rows = 0;
        maxCols = 0;
        buf.clear();
        count.clear();
        data = 0;
    }

    const IndexedValue* rowPtr(int row) const
    {
        CV_Assert(row >= 0 && row < rows);
        return data + row * maxCols;
    }

    IndexedValue* rowPtr(int row)
    {
        CV_Assert(row >= 0 && row < rows);
        return data + row * maxCols;
    }

    void insert(int row, int col, double value)
    {
        CV_Assert(row >= 0 && row < rows);

        int currCount = count[row];
        CV_Assert(currCount < maxCols);

        IndexedValue* rowData = rowPtr(row)