本次实验的基础是在《Fast and Effective L0 Gradient Minimization by Region Fusion》和《Feature-preserving filtering with L0 gradient minimization 》两篇论文的基础上进行的改进。
这两篇论文中都是利用合并平均值相同的区域的思想来实现L0范式最小化以及保持图像的边缘特性。在这两篇论文中都用到了将L0范式转化为L2范式的思想。
本次实验的基本思想分为两种：
第一种：利用论文种提到的算法，直接将原始图片进行流场平滑，在对平滑后的图片进行流场绘制。
第二种：利用论文中提到的算法，首先将原始图片转化为流场，再将流场看作一幅待处理的原始图片进行平滑。
方案一的实现效果
原始图片
参数为0.1
参数为0.2
参数为0.3
参数为0.6
参数为1
参数为1.5
参数为1.7
参数为1.8
方案二的实现效果
原始图片
参数为0.1
参数为0.2
参数为0.3
参数为0.4
参数为0.5
参数为20
代码
将图片转化为流场
Mat image_to_flow_field(Mat & image)
{
//创建grad_x和grad_y矩阵
Mat grad_x = Mat::zeros(image.size(), CV_16S);
Mat grad_y = Mat::zeros(image.size(), CV_16S);
//创建roberts算子处理显示图片
Mat roberts_img = Mat::zeros(image.size(), CV_16S);
for (int i = 0; i < image.rows - 1; i++)
{
for (int j = 0; j < image.cols - 1; j++)
{
//求取roberts算子水平方向的梯度
grad_x.at<short>(i, j) = image.at<uchar>(i, j) - image.at<uchar>(i + 1, j + 1);
//求取roberts算子垂直方向的梯度
grad_y.at<short>(i, j) = image.at<uchar>(i, j + 1) - image.at<uchar>(i + 1, j);
}
}
//构建流场
for (int i = 0; i < image.rows; i++)
{
for (int j = 0; j < image.cols; j++)
{
if (grad_x.at<ushort>(i, j) != 0)
{
//求每个像素的角度
roberts_img.at<short>(i, j) = atan(grad_y.at<short>(i, j) / grad_x.at<short>(i, j)) * 180 / 3.1415;
}
else
{
//计算特殊角度
roberts_img.at<short>(i, j) = 90;
}
}
}
roberts_img.convertTo(roberts_img, CV_32FC1);
return roberts_img;
}L0测度的流场平滑
Mat l0_norm(Mat img, double lambda, int maxSize, int maxLoop)
{
//获取图像的基本属性：长、宽、通道数
int rows = img.rows, cols = img.cols, bands = img.channels();
//获取图像的大小
int M = rows * cols;
//定义基本变量
double total = 0;
clock_t start, end;
double time;
double step = (double)1/maxLoop;
double ct = 0;
double curThresh;
//定义存储图片各个组中元素均值的数组
float* Y = new float[bands*M];
//定义存储图片各个组中元素和的数组
float* sumY = new float[bands*M];
//将原始图片的数据赋值给Y和sumY
for(int i = 0; i<bands*M; i++)
{
sumY[i] = Y[i] = img.data[i];
}
//定义数量为M（原始图片像素个数）的组链表
LinkedList* G = new LinkedList[M];
//定义存储每个元素相邻四个元素的数组
int* NB = new int[M*maxSize];
//定义存储每个元素中相邻元素的数目，注意存储的是实际数目的2倍
int* nNB = new int[M];
//定义每组中的最大元素值
int* maxNB = new int[M];
//定义每组中元素的个数
int* W = new int[M];
//定义每组的索引值
int* IDX = new int[M];
int* lutIDX = new int[M];
//定义每组的hash索引
int* hashID = new int[M];
//定义每组的默认页
bool* pagefault = new bool[M];
//定义每组的起始页
int** startpage = new int*[M];
for(int i = 0; i < M; i++)
{
//构建每个组的链表
G[i].insert(i);
//构建每个组的索引
IDX[i] = i;
lutIDX[i] = i;
//构建每个组的初始数量为1
W[i] = 1;
//构建每个组的最大容量
maxNB[i] = maxSize;
//构建每个组的hash索引
hashID[i] = -1;
//构建每个组的页的状态
pagefault[i] = false;
//构建每个页的起始地址
startpage[i] = &NB[i*maxSize];
}
//创建每个组的邻居节点
createNeighbour(rows, cols, maxSize, NB, nNB);
//定义迭代变量
int iter = 0;
//循环迭代
while(iter <= maxLoop && M > 1)
{
int maxNBnum = 0;
//double curThresh = pow(1.5,iter-maxLoop)*t;
//double curThresh = pow(ct,2.2)*t;
//定义辅助变量的线性增长函数
double curThresh = ct * lambda;
//辅助变量beta逐步增加
ct += step;
for(int i = 0; i < M; i++)
{
int idx1 = IDX[i];
int* startI1 = startpage[idx1];
//获取每组中的元素的数目
int FIX_LOOP_TIMES = nNB[idx1];
//创建各个元素的hash索引
for(int hi = 0; hi < nNB[idx1]; hi = hi+2)
{
hashID[*(startI1+hi)] = hi;
}
//依次遍历每个组的元素
for(int j = 0; j < FIX_LOOP_TIMES; j=j+2)
{
//获取每组的起始地址
int idx2 = *(startI1+j);
int* startI2 = startpage[idx2];
int rIdx1 = bands*idx1;
int rIdx2 = bands*idx2;
int len = *(startI1+j+1);
//计算||Yi-Yj||2
float dx = Y[rIdx1 ] - Y[rIdx2 ];
float dy = Y[rIdx1+1] - Y[rIdx2+1];
float dz = Y[rIdx1+2] - Y[rIdx2+2];
double d = dx*dx + dy*dy + dz*dz;
//计算相邻两个组的元素个数之和
int sumW = W[idx1]+W[idx2];
//计算公式wiwj||Yi −Yj||2 ≤ βci,j(wi + wj)
if(d*W[idx1]*W[idx2] <= curThresh*len*sumW)
{
//合并相邻的两个组，并计算均值
for(int b = 0; b < bands; b++)
{
//将索引为2的组中的所有元素和累加到索引为1的组中
sumY[rIdx1+b] += sumY[rIdx2+b];
//重新求取索引为1的组的元素均值
Y[rIdx1+b] = sumY[rIdx1+b]/sumW;
}
//更新合并后的组的元素数目以及组的元素
W[idx1] = sumW;
G[idx1].append(G[idx2]);
//从idx1中删除idx2
if(j != nNB[idx1]-2)
{
swap(*(startI1+j), *(startI1+nNB[idx1]-2));
swap(*(startI1+j+1), *(startI1+nNB[idx1]-1));
hashID[*(startI1+j)] = j;
j = j - 2;
}
nNB[idx1] -= 2;
//移除一个元素后总次数减少2个。
if(nNB[idx1] < FIX_LOOP_TIMES)
FIX_LOOP_TIMES -= 2;
//更新hash表
hashID[idx2] = -1;
//将idx2的相邻节点插入到idx1
for(int t =0 ; t<nNB[idx2]; t = t+2)
{
int aa = *(startI2+t);
int la = *(startI2+t+1);
int* startIa = startpage[aa];
if (aa == idx1)
continue;
//查看aa是否有idx1和idx2的相同邻域
int find = hashID[aa];
//如果idx1和idx2有相同的邻域
if (find > -1)
{
*(startI1+find+1) += la;
int k = 0;
while(*(startIa+k) != idx1)
{
k += 2;
}
//更新新的连接数目
*(startIa+k+1) = *(startI1+find+1);
//从idx2的邻域中清除idx2
k = 0;
while(*(startIa+k) != idx2)
{
k += 2;
}
if(k != nNB[aa]-2)
{
swap(*(startIa+k ), *(startIa+nNB[aa]-2));
swap(*(startIa+k+1), *(startIa+nNB[aa]-1));
}
nNB[aa] -= 2;
}
如果idx1和idx2没有相同的邻域
else
{
if(nNB[idx1] >= maxNB[idx1])
{
if(pagefault[idx1] || maxNB[idx1+maxNB[idx1]/maxSize] != 0)
{
// PAGE FAULT
int*temp = new int [2*maxNB[idx1]];
maxNB[idx1] = 2*maxNB[idx1];
// Copy to new page
for(int ii = 0; ii < nNB[idx1]; ii++)
{
temp[ii] = *(startI1+ii);
}
if(pagefault[idx1])
{
delete[] startpage[idx1];
}
startpage[idx1] = &temp[0];
startI1 = startpage[idx1];
pagefault[idx1] = true;
}
else
{
maxNB[idx1] += maxSize;
}
}
//更改C(i,j)
*(startI1+nNB[idx1]) = aa;
*(startI1+nNB[idx1]+1) = la;
//更新hash表
hashID[aa] = nNB[idx1];
nNB[idx1] += 2;
int k = 0;
while(*(startIa+k) != idx2)
{
k += 2;
}
*(startIa+k) = idx1;
}
}
// DELETE!
maxNB[idx2] = 0;
nNB[idx2] = 0;
if(pagefault[idx2])
delete[] startpage[idx2];
M = M - 1;
int p_idx2 = lutIDX[idx2];
lutIDX[IDX[M]] = p_idx2;
swap(IDX[p_idx2], IDX[M]);
}
}
//更新idx1的最大值
if (nNB[idx1] > maxNBnum)
maxNBnum = nNB[idx1];
//清除hashIDX
for(int hi = 0; hi < nNB[idx1]; hi = hi+2)
{
hashID[*(startI1+hi)] = -1;
}
}
//输出原始图像中的组数
if (iter == 0)
{
cout << "Original number of groups = " << M << endl;
}
//迭代次数加1
iter++;
}
//输出平滑后的最终组数目
cout<<"Final number of groups = "<<M<<endl;
// restore image
for(int i = 0; i < M; i++)
{
int idx1 = IDX[i];
if(pagefault[idx1])
{
delete[] startpage[idx1];
}
Node* temp = G[idx1].pHead;
while(temp != NULL)
{
int lidx = bands*temp->value;
int ridx = bands*idx1;
for(int b = 0; b < bands; b++)
img.data[lidx+b] = Y[ridx+b];
temp = temp->next;
}
}
//清空所占据的内存
delete[] Y;
delete [] NB;
delete [] nNB;
delete[] W;
delete[] sumY;
delete[] IDX;
delete[] lutIDX;
delete[] hashID;
delete[] G;
delete[] pagefault;
delete[] startpage;
//返回平滑后的图像
return img;
}总结
- 通过此次实验我发现针对流场的平滑并不能直接利用现有的、较为成熟的算法来应用于流场的平滑。
- 针对两种不同的思路：一种是直接平滑原始图片，一种是直接平滑流场。通过对比发现，利用直接平滑图片，间接平滑流场这种方式，取得的效果较好。但是，这实际上就是将问题转化为了平滑图像，这与我实际的出发点并不相符合，因此，接下来我打算直接利用平滑流场的方式来进行处理。
- 在实验的过程中，我发现大多数论文中仅仅提及到了一些其算法的基本思想，如果多看几遍可以了解他的思想，但是但想要实现作者提出的算法时，就会发现算法中有很多的细节以及关键参数的选取，在论文中并没有提及。这时，你即使联系这篇论文的第一作者或通讯作者，也基本上得不到任何回复。我经过亲自尝试无果后，发现可以根据每篇论文的参考文献来查找相关的线索，一般情况下，要想真正的理解一篇论文的核心思想，就至少要看10篇相关的论文。换句话说，只要把一篇论文的参考文献都看完就可以了。
- 我发现大多数论文的算法根本就不具备通用性，这些算法可能仅仅针对自己研究领域的某一小类图像处理效果较好，但是针对其他图像的处理效果较差。
- 另外，几乎所有的论文都是不公开源代码的，这就无法做到实验再现！这不禁让我怀疑这个实验效果是不是真的。在国外的一些论文中，在很多大牛的个人主页都公布了相关的论文代码，我认为这是值得学习的地方。
- 在前段时间，我花费了很长的时间来查找，下载文献，主要是因为有很多论文都需要花钱买。在这里我给大家推荐一个比较好的文献集：“80图书馆”。这里包括了绝大部分英文、中文的常用数据库。
467
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
