视差图反向合成虚拟视点/backward view synthesis base on disparity

视差图的反向合成算法，就是根据已知的虚拟视点的视差图和对应的已知视点的原图进行合成虚拟视点的过程。假设使用middlebury数据集为例，已知im2和对应的视差图disp6，我们现在要求im6的图像，就可以根据视差图反向合成虚拟视点这一方法进行求取。反向合成过程，我们已知一点的水平坐标假设为1，对应的虚拟视点的视差为4，那么我们根据公式 Xl - Xr = d 可以知道的是 Xl = d + Xr为5，所以我们在已知视图的位置进行取对应的像素。但是这个过程中可能会出现Xl不是整数的情况，比如计算的Xl是4.5之类的，这个时候可以采用双线性插值或者最近邻插值的方法。

还有一种就是根据已知的视图和对应的视差图，先将视差图进行前向合成得到虚拟视点的视差图，再根据得到的虚拟视点的视差图在已知的视图中进行取值得到虚拟视点的图像。但是在这过程中，因为前向合成得到的视差图存在破洞，所以要对破洞进行处理。

下面是对应的第二种方法进行虚拟视点合成的算法：

#include <iostream>
#include <string>
#include <opencv.hpp>
#include <vector>
#include <opencv2/highgui/highgui_c.h>  

using namespace std;
using namespace cv;

int index(int m, int n)
{   // 如果没有超过边界 直接返回对应的值 
    if (m>=0 && m<n)
        return m;
    // 如过超过了左边界 那么直接返回0 我觉得不太合适 下面的也是 具体再看看
    else if (m<0)
        return 0;
    // 如果超过了右边界 直接返回 最右边的列 
    else if (m>=n)
        return n-1;
    return 0;
}

/* 
 天降大任于斯人也，必先劳其筋骨，饿其体肤
*/

void obtainNewDispMap(const Mat &refDisp, Mat &dstDisp, float value)
{
    int height=refDisp.rows;
    int width=refDisp.cols;
    uchar* pSrcDispData=(uchar*) refDisp.data;
    float* pDstDispData=(float*) dstDisp.data;
    for (int j=0; j<height; j++)
    {
        for (int i=0; i<width; i++)
        {
            int disp=(int)pSrcDispData[j*width+i]; // 取源视差图的视差值 并强转为 int
            float newDisp=disp*(value);  // 这value 为1  一般设为 
            int inew=(int)(i-newDisp);  // 根据视差的值 获取对应 目标视差图的位置 
            inew=index(inew, width);  // 判断这个位置是否超出了范围 
            if(uchar(newDisp) >= 255)
                newDisp = 0;
            pDstDispData[j*width+inew]=newDisp;  // 根据获取的目标视差图的位置 inew 设置 那一个坐标的视差为 newDisp 
        }
    }
}


void left_and_right_consistance(Mat left_disp, Mat right_disp, vector<std::pair<int, int>> occlusions_, 
    vector<std::pair<int, int>> mismatch_)
{
    occlusions_.clear();
    mismatch_.clear();
    int row = left_disp.rows;
    int col = left_disp.cols;
    for(int i=0; i<row; i++)
    {
        for(int j=0; j<col; j++)
        {
            auto disp_l = left_disp.at<uchar>(i, j);
            auto col_right = uchar(j - disp_l + 0.5);
            if(col_right >=0 && col_right < col)
            {
                auto disp_r = right_disp.at<uchar>(i, col_right);
                auto col_left = uchar(col_right + disp_r + 0.5);
                if(col_left >=0 && col_left < col)
                {
                    auto disp_l_ = left_disp.at<uchar>(i, col_left);
                    if(disp_l_ > disp_l)
                        occlusions_.emplace_back(i, j);
                    else
                        mismatch_.emplace_back(i, j);
                }
            }
        }
    }
}

/*
* param:disp 需要填充的mat 
* return: void 
* 直接是图上左边的值 传给空的位置 
*/
void fillDispHole(Mat& disp)
{
    int col = disp.cols;
    int row = disp.rows;
    for(int i=0; i<row; i++)
    {
        for(int j=0; j<col; j++)
        {
            float_t disp_val = disp.at<float_t>(i, j);
            if(disp_val == 0 || disp_val == 255)
            {
                int index = j-1;
                while(index>=0 && (disp.at<float_t>(i, index) == 0 || disp.at<float_t>(i, index) == 0))
                {
                    index--;
                }
                if(index==0)
                    continue;
                disp.at<float_t>(i, j) = disp.at<float_t>(i, index);
            }
        }
    }
}


/*
The most teribble thing is cancer
I need prepare for it 
first thing is money
second thing is my paper include experiments
third is accompany, care and comfort
*/

// 5000 + 5000  
// 7400 + 15000 + 3000 = 25400 

/*
* param: disp 需要调整对比度的视差图
* return：void 因为是引用传值 
*/
void AdjustContrast(Mat disp)
{
    double minval = 0.0;
    double maxVal = 0.0;
    minMaxLoc(disp, &minval, &maxVal);
    for(int i=0; i<disp.rows; i++)
    {
        for(int j=0; j<disp.cols; j++)
        {
            uchar val = disp.at<uchar>(i, j);
            double adjust_val = (double(val) - minval) / (maxVal - minval) * 255;
            // cout << uchar(adjust_val) << endl;
            disp.at<uchar>(i, j) = uchar(adjust_val);
        }
    }
}


void insertDepth32f(cv::Mat& depth)
{
    const int width = depth.cols;
    const int height = depth.rows;
    float* data = (float*)depth.data;
    cv::Mat integralMap = cv::Mat::zeros(height, width, CV_64F);
    cv::Mat ptsMap = cv::Mat::zeros(height, width, CV_32S);
    double* integral = (double*)integralMap.data;
    int* ptsIntegral = (int*)ptsMap.data;
    memset(integral, 0, sizeof(double) * width * height);
    memset(ptsIntegral, 0, sizeof(int) * width * height);
    for (int i = 0; i < height; ++i)
    {
        int id1 = i * width;
        for (int j = 0; j < width; ++j)
        {
            int id2 = id1 + j;
            if (data[id2] > 1e-3)
            {
                integral[id2] = data[id2];
                ptsIntegral[id2] = 1;
            }
        }
    }
    // 积分区间
    for (int i = 0; i < height; ++i)
    {
        int id1 = i * width;
        for (int j = 1; j < width; ++j)
        {
            int id2 = id1 + j;
            integral[id2] += integral[id2 - 1];
            ptsIntegral[id2] += ptsIntegral[id2 - 1];
        }
    }
    for (int i = 1; i < height; ++i)
    {
        int id1 = i * width;
        for (int j = 0; j < width; ++j)
        {
            int id2 = id1 + j;
            integral[id2] += integral[id2 - width];
            ptsIntegral[id2] += ptsIntegral[id2 - width];
        }
    }
    int wnd;
    double dWnd = 2;
    while (dWnd > 1)
    {
        wnd = int(dWnd);
        dWnd /= 2;
        for (int i = 0; i < height; ++i)
        {
            int id1 = i * width;
            for (int j = 0; j < width; ++j)
            {
                int id2 = id1 + j;
                int left = j - wnd - 1;
                int right = j + wnd;
                int top = i - wnd - 1;
                int bot = i + wnd;
                left = max(0, left);
                right = min(right, width - 1);
                top = max(0, top);
                bot = min(bot, height - 1);
                int dx = right - left;
                int dy = (bot - top) * width;
                int idLeftTop = top * width + left;
                int idRightTop = idLeftTop + dx;
                int idLeftBot = idLeftTop + dy;
                int idRightBot = idLeftBot + dx;
                int ptsCnt = ptsIntegral[idRightBot] + ptsIntegral[idLeftTop] - (ptsIntegral[idLeftBot] + ptsIntegral[idRightTop]);
                double sumGray = integral[idRightBot] + integral[idLeftTop] - (integral[idLeftBot] + integral[idRightTop]);
                if (ptsCnt <= 0)
                {
                    continue;
                }
                data[id2] = float(sumGray / ptsCnt);
            }
        }
        int s = wnd / 2 * 2 + 1;
        if (s > 201)
        {
            s = 201;
        }
        cv::GaussianBlur(depth, depth, cv::Size(s, s), s, s);
    }
}


int main()
{
    // 读数据 
    string imgPath="middlebury\\teddy\\";
    Mat srcImgL=imread(imgPath+"im2.png");
    Mat dispL=imread(imgPath+"disp2.png", 0);
    dispL=dispL/4;  // 因为是1/4分辨率 

    int imgHeight=srcImgL.rows;
    int imgWidth=srcImgL.cols;

    // 设置容器 存储
    Mat dstImgL=Mat::zeros(imgHeight,imgWidth, CV_8UC3);
    Mat dstImg=Mat::zeros(imgHeight,imgWidth, CV_8UC3);
    Mat dstNewDispImg=Mat::zeros(imgHeight,imgWidth, CV_32FC1);
    Mat save_disp = Mat::zeros(imgHeight, imgWidth, CV_8UC1);

    // 初始化指针变量
    uchar* pImgDataL=(uchar*)srcImgL.data;
    uchar* pDispDataL=(uchar*)dispL.data;
    uchar* pDstDataL=(uchar*)dstImgL.data;

    dstNewDispImg.setTo(0);
    dstImgL.setTo(0);
    save_disp.setTo(0);
    // if (cnt%2==0)   interp=(float)k/viewCnt;
    // else interp=(float)(viewCnt-k)/viewCnt;

    // 获取虚拟视点的视差图 dstNewDispImg
    obtainNewDispMap(dispL, dstNewDispImg, 1);
    insertDepth32f(dstNewDispImg);
    // fillDispHole(dstNewDispImg);

    // double minVal1 = 0.0;
    // double maxVal1 = 0.0;
    // minMaxLoc(dstNewDispImg, &minVal1, &maxVal1);
    // int maxval_int = int(maxVal1);
    // int minval_int = int(minVal1);
    // cout << "minVal1:" << minval_int << endl;
    // cout << "maxVal1:" << maxval_int << endl;

    float* pNewDispData=(float*)dstNewDispImg.data;
    for (int j=0; j<imgHeight; j++)
    {
        for (int i=0; i<imgWidth; i++)
        {   
            // 获取新视差图的视差
            float disp=pNewDispData[j*imgWidth+i];
            float id=i+disp;  // 获取目标视图的位置 


            int id0=floor(id);  // 向下取值  加入坐标是 1.7 那么获得id0 就是 1 id1 就是 2
            int id1=floor(id+1);


            //if(id0 == 255 || id1 == 255)
            //    cout << " id = " << id << endl;

            float weight1=1-(id-id0);  // weight 是 id0的权重 同上 那么weight1的权重为 0.3 weight2为 0.7
            float weight2=id-id0;

            // 判断是否超出了边界
            id0=index(id0, imgWidth); 
            id1=index(id1, imgWidth);


            save_disp.at<uchar>(j, i) = uchar(disp);

            //插值结果
            pDstDataL[j*imgWidth*3+i*3+0]=weight1*pImgDataL[j*imgWidth*3+id0*3+0]+weight2*pImgDataL[j*imgWidth*3+id1*3+0];
            pDstDataL[j*imgWidth*3+i*3+1]=weight1*pImgDataL[j*imgWidth*3+id0*3+1]+weight2*pImgDataL[j*imgWidth*3+id1*3+1];
            pDstDataL[j*imgWidth*3+i*3+2]=weight1*pImgDataL[j*imgWidth*3+id0*3+2]+weight2*pImgDataL[j*imgWidth*3+id1*3+2];
        }
    }
    // fillDispHole(save_disp);
    // namedWindow("virImg");
    // imshow("virImg", dstImgL);
    // dstNewDispImg = int(dstNewDispImg);
    AdjustContrast(save_disp);
    imwrite("save_disp.png", save_disp);
    imwrite("img_syn_backward.png", dstImgL);
    // waitKey(10);
    return 0;
}

输入的图像是middlebury的teddy数据集：如下所示分别是im2和对应的disp2：

 

得到的结果如下所示;

 

如结果所示，得到的图像不存在空洞问题。这是反向合成相比前向合成的优势。

但是因为对不存在的位置进行了前向求取视点，所以导致右下角的部分其实不是真实的，所以对应产生了伪影的问题。

middlebury数据集的网址如下所示：

2003 Stereo Datasets