插图.jpg
一. 前言：
当年在武大电信院做毕业设计的时候，我的导师希望我用C代码完成彩色图像的雾霾程度评价，如果时间允许，能够开发出一款APP完成移动端的实现。很可惜，当时我对于iOS开发一无所知，时间不允许我完成移动端的实现，但是完成了Windows端的实现，利用MFC做了一个毕业设计展示。还记得当时花了很多时间研究图像处理，一步一步实现每一个算法细节，辛苦却可以感受到满满的收获。
第一份工作找的并不是开发，是软件测试，不过这也给了我时间去完善我当年毕设没有实现移动端遗憾。花了很长时间利用下班时间自学iOS,技术还是很渣，希望大家能够勉强看得下去我的文章。
图像中雾霾程度的评价是图像增强处理的第一步，从方法上可以分为主观评价方法和客观评价方法，前者凭感知者主观感受评价;后者依据模型给出的量化指标或参数衡量。客观评价的目标使得评价模型准确地反映人眼视觉感知的主观评价。本课题研究基于自然场景统计特征(NSS)和其他感知特征的雾霾浓度评价研究，使得对于输入的图像，准确输出评价值。
先给展示一下demo
demo.gif
操作步骤：
1.选取图片，然后可以从手机相册中选择一张图片，也可以拍照
2.可以修改尺寸大小，点击完成
3.点击雾霾程度分析的button，就可以看到给出这张图片的雾霾程度的客观评价值
二. 前期准备
考虑到 OpenCV 是基于 C/C++ 可跨平台的通用 Lib，为了降低学习成本，便将整个学习和实验集成到 iOS 的开发环境里了。前期要做如下几方面的准备工作：
直接到OpenCV官网下载你想要的OpenCV的framework；
将 OpenCV.framework 导入 iOS 项目工程中；
因为 OpenCV 中的 MIN 宏和 UIKit 的 MIN 宏有冲突，所以需要在 .pch 文件中，先定义 OpenCV 的头文件，否则会有编译错误；
将需要混编 C++ 和 Objective-C 的文件后缀改为 .mm;
为 UIImage 添加 Category，方便与OpenCV 图象格式的数据 cv::Mat 相互转换。
因这些繁琐的配置问题不是本文写作重点，而且网上不乏一些详细说明，推荐参考唐巧大神的 在MacOS和iOS系统中使用OpenCV 一文，这里就不再赘述。
三. 算法：基于自然场景统计特征（NSS）和其他的雾感统计特征的雾霾程度分析
下面是算法的实现流程图：
image
简单介绍下个这个算法：
空域 NSS 模型即计算自然图像上的局部差值均值和对比度归一化(MSCN) 系数。雾感统计特征包括 MSCN 系数，局部平均和局部系数的锐利度方差，对比能，图像熵，像素暗通道先验，色彩饱和度和色彩丰富度等，总共 12个雾感知统计特征，利用这些雾感统计特征去计算每一个已经被分割为 P×P 大小 的图像块，然后得到图像的多元高斯(MVG)模型，利用 Mahalanobis 类距离测量方法，测量图像的雾感知统计特征的 MVG 模型分别与自然雾霾图像库的 MVG模型和自然无雾图像库的 MVG 模型之间的距离，最后计算得到感知雾霾程度。
三.算法具体实现细节和主要代码
这里仅贴一下部分我认为比较重要的代码，具体实现细节请参考我的github：https://github.com/echooj/FogDetect
当然，懒得看代码的同学可以直接跳到测试结果部分，就是有一张可爱的小蘑菇图的位置
基础过程
图像尺寸的测量和 RGB 值的获取
先提取图像的尺寸，然后修改图像尺寸，使得图像的宽和高能够被8整除，最后提取修改后的图像的R、G、B值和灰度值Ig。
vector<Mat> channels; Mat Ig; Mat Blue; Mat Green; Mat Red; cvtColor(I2, Ig, CV_BGR2GRAY); split(I2, channels); Blue = channels.at(0); Green = channels.at(1); Red = channels.at(2);
HSV 空间的转换
将图像转换到HSV空间，并获取s(饱和度)的值。
Mat I_hsv; Mat Is;vector<Mat> channels2; cvtColor(I2, I_hsv, CV_BGR2HSV); split(I_hsv, channels2); Is = channels2.at(1); Is.convertTo(Is,CV_64F); Is=Is/255;
暗通道先验
利用了何凯明的暗通道先验的实现过程。首先先将R,G,B格式转换从uchar到double，然后进行 归一化处理，使得取值范围在0-1之间，最后利用两个min得到Id。
Mat blue = Mat_<double>(Blue); Mat green = Mat_<double>(Green); Mat red = Mat_<double>(Red); Mat Irn = red / 255; Mat Ign = green / 255; Mat Ibn = blue / 255; Mat Id; min(Irn, Ign, Id); min(Id, Ibn, Id);
MSCN 系数的获得以及 rg 和 by 通道的获得
获得MSCN系数和cv。先制作MSCN_window的模板，然后对Ig和 MSCN_window进行实现线性空间滤波函数中的相关运算得到mu;同样利用Ig.*Ig和 MSCN_window的滤波运算再和mu_sq进行做差，最后开方得到sigma;利用(Ig- mu)./(sigma+1)得到MSCN系数。rg和by通道获得是利用OpenCV函数进行常规的矩 阵运算即可。
// MSCNMat Ignew; Ig.convertTo(Ignew, CV_64F); Mat mu = Mat::zeros(Ig.rows, Ig.cols, CV_64FC1); mu = conv2(Ignew, MSCN_window, CONVOLUTION_SAME); Mat DIg = Ignew.mul(Ignew); Mat Mattemp = conv2(DIg, MSCN_window, CONVOLUTION_SAME); Mat sigma; Mat Mattemp1; Mat mu_sq = Mat::zeros(mu.rows, mu.cols, CV_64FC1); mu_sq = mu.mul(mu); absdiff(Mattemp, mu_sq, Mattemp1); sqrt(Mattemp1, sigma); Mat temp1 = Ignew - mu; Mat temp2 = sigma + 1; Mat MSCN; MSCN = temp1 / temp2; Mat cv = sigma / mu;// rg 和 by 通道Mat R, G, B; Red.convertTo(R, CV_64F); Green.convertTo(G, CV_64F); Blue.convertTo(B, CV_64F); Mat rg = R - G; Mat by = 0.5*(R + G) - B;
雾感知统计特征的提取
f1,f2...f12共12个雾感统计特征，具体功能可以参考表格2.1
image.png
f1:
(1)概述:f1的提取，f1的功能可以参考表格2.1。
Mat tempbb = im2col(MSCN, ps); Mat tempcc = var(tempbb);//Mat MSCN_var = NewReshape(tempcc, row / ps);
(2)im2col的功能:重排矩阵列;
Mat im2col(Mat &InputMat, int ps){ int row = InputMat.rows / ps;int col = InputMat.cols / ps;
Mat **t = new Mat *[row];for (int i = 0; i<row; ++i)
t[i] = new Mat[col];
Mat **t1 = new Mat *[row];for (int i = 0; i<row; ++i)
t1[i] = new Mat[col];
Mat **xin = new Mat *[row];for (int i = 0; i<row; ++i)
xin[i] = new Mat[col];for (int i = 0; i < col; I++)for (int j = 0; j < row; j++)
{ xin[j][i] = InputMat(Rect(i*ps, j*ps, ps, ps));
t[j][i] = xin[j][i].clone();
Mat temp = t[j][i].t();
t1[j][i] = temp.reshape(0, ps*ps);}
Mat outputMat = Mat::zeros(ps*ps, row*col, CV_64F);
Mat *t2 = new Mat[row*col];int i, j, count = 0;for (i = 0; i < col; I++)
{for (j = 0; j < row; j++)
{
t2[count] = t1[j][I];
count++;
}
}for (int i = 0; i < outputMat.rows; I++)for (int j = 0; j < outputMat.cols; j++)
outputMat.at<double>(i, j) = t2[j].at<double>(i, 0);return outputMat;}(3)nanvar的功能:忽视掉NaN的方差的计算，注意分母是n-1;
(4)reshape的功能:重新调整矩阵的行数、列数。
Mat NewReshape(Mat &InputMat, int ps){ int row = 1;int col = InputMat.cols / ps;
Mat *t = new Mat[col];
Mat *t1 = new Mat[col];
Mat *xin = new Mat[col];for (int i = 0; i < col; I++)
{ xin[i] = InputMat(Rect(i*ps, 0, ps, 1));
t[i] = xin[i].clone();
t1[i] = t[i].reshape(0, ps);}
Mat outputMat = Mat::zeros(ps, row*col, CV_64F);for (int i = 0; i < outputMat.rows; I++)for (int j = 0; j < outputMat.cols; j++)
outputMat.at<double>(i, j) = t1[j].at<double>(i, 0);return outputMat;}雾感知统计特征中 f2 和 f3 的提取
(1)概述:f2和f3的提取，详情请参考以下代码。
Mat temp_vertical = circshift(MSCN, 1); Mat temp_vertical2 = MSCN.mul(temp_vertical); Mat MSCN_V_pair_col = im2col(temp_vertical2, ps); Mat MSCN_V_pair_col_temp1 = MSCN_V_pair_col.clone();for (int i = 0; i < MSCN_V_pair_col_temp1.rows; I++)for (int j = 0; j < MSCN_V_pair_col_temp1.cols; j++)if (MSCN_V_pair_col_temp1.at<double>(i, j)>0) MSCN_V_pair_col_temp1.at<double>(i, j) = 0; Mat MSCN_V_pair_col_temp2 = MSCN_V_pair_col.clone();for (int i = 0; i < MSCN_V_pair_col_temp2.rows; I++)for (int j = 0; j < MSCN_V_pair_col_temp2.cols; j++)if (MSCN_V_pair_col_temp2.at<double>(i, j)<0) MSCN_V_pair_col_temp2.at<double>(i, j) = 0; Mat MSCN_V_pair_L_var=NewReshape(nanvar(MSCN_V_pair_col_temp1),row/ps); Mat MSCN_V_pair_R_var=NewReshape(nanvar(MSCN_V_pair_col_temp2),row/ ps);
(2)cirshift: 矩阵循环平移。因为整个项目只用到了一次循环移位，所以该函数是特殊的 循环移位，也就是按列方向循环向下移动一位的函数。
雾感知统计特征的 f4 和 f5 的提取
(1)概述:f4和f5的提取，使用的算法和f2、f3相似。即reshape配合mean和 im2col的函数。
(2)mean:获取矩阵中每一列的平均值。雾感知统计特征的 f6、f7 和 f8 的提取
(1)概述:本模块包括CE的编写，CE中又使用了border_in和border_out函数。
(2)CE:感知对比能，获取灰度、蓝黄、红绿颜色通道。
(3)border_in: 边界增加图像。
(4)border_out:边界平衡图像。雾感知统计特征的 f9 、f10 和 f11 的提取
(1)概述:f9、f10、f11的提取。
(2)Entropy:实现图像熵的功能。
Entropy的即图像熵-sum(p.*log2(p))的实现算法:
a. 利用temp[256]数组做一个类似统计表的功能，统计图像对应的0-255灰度值对 应的像素个数;
b. 计算每一个个像素的概率;
c. 根据定义计算图像熵。
double Entropy(Mat img){ double temp[256];for (int i = 0; i<256; I++)
temp[i] = 0.0;// 计算每个像素的累积值for (int m = 0; m<img.rows; m++)
{ double *t = img.ptr<double>(m);for (int n = 0; n<img.cols; n++)
{ int i = t[n];
temp[i] = temp[i] + 1;}
}// 计算每个像素的概率for (int i = 0; i<256; I++)
{temp[i] = temp[i] / (img.rows * img.cols);}double result = 0;// 计算图像熵for (int i = 0; i<256; I++)
{ if (temp[i] == 0.0)
result = result;elseresult = result - temp[i] * (log(temp[i]) / log(2.0));
}return result; }(3)num2cell:将矩阵中的数据转成一个一个的孢元。
(4)cellfun:对每个孢元单独计算.
雾感知统计特征的 f12 的提取
功能
概述:对以上f1—f12所有的特征进行重排矩阵列，然后放到一个矩阵feat中, 最后对feat进行变化得所要的值。
MVG 模型的距离计算
概述:MVG模型Df的计算和Dff的计算，下面将主要介绍Df的计算，Dff的计算 除了数据不同外，基本上就是Df的重复。
自然无雾图像特征的提取和块的测试参数提取
自然无雾图像特征的提取:载入两个矩阵，cov_fogfreeparam和 mu_fogfreeparam块的测试参数提取:对之前的feat进行一些信息的提取。距离计算
(1)概述:实现距离计算，包括了很多中间步骤雾霾等级计算
(1)概述:Df和Df_map的获得，这个部分的完结。
Mat** distance_patch_t1 = new Mat *[mu_matrix.rows];for (int i = 0; i<mu_matrix.rows; ++i)
distance_patch_t1[i] = new Mat[1];for (int i = 0; i < mu_matrix.rows; I++)for (int j = 0; j < 1; j++)
distance_patch_t1[i][j] = mu_cell[i] * (cov_cell[i][j].inv());
Mat** distance_patch_t2 = new Mat *[mu_matrix.rows];for (int i = 0; i<mu_matrix.rows; ++i)
distance_patch_t2[i] = new Mat[1];for (int i = 0; i < mu_matrix.rows; I++)for (int j = 0; j < 1; j++)
distance_patch_t2[i][j]=distance_patch_t1[i][j]* (mu_transpose_cell[I]);
Mat distance_patch =
Mat::zeros(distance_patch_t2[0][0].rows*mu_matrix.rows,
distance_patch_t2[0][0].cols, CV_64F);for (int i = 0; i < distance_patch.rows; I++)for (int j = 0; j < distance_patch.cols; j++)
{
distance_patch.at<double>(i,j)=distance_patch_t2[i][0].at<double>(0,j);
}sqrt(distance_patch, distance_patch);
Mat Df = NewMean(distance_patch);
Mat distance_patch_t = distance_patch.t();
Mat Df_map = NewReshape(distance_patch_t, row / ps);雾霾程度的获取
概述:首先获取雾霾密度值，然后根据雾霾密度值计算雾霾程度。即先利用D = Df / (Dff + 1)获取雾霾密度值，然后根据雾霾密度值D在0 ≤ D < 1为低，1 < D ≤ 3为中，雾霾密度值D在D > 3为高，进行雾霾程度输出。
四. 测试结果
首先我们载入一张雾霾程度比较低的图片
作者：linlin泠迹
链接：https://www.jianshu.com/p/e44c26b7eb3e
共同学习，写下你的评论
评论加载中...
作者其他优质文章
