图像增强--视网膜皮层Retinex算法（二）

视网膜皮层Retinex算法可以参考下面两篇文章

http://www.cnblogs.com/Imageshop/archive/2013/04/17/3026881.html

https://blog.csdn.net/carson2005/article/details/9502053

下面直接进行讲解基于Opencv 4.10的Retinex算法实现

/*主函数*/
int main()
{Mat orig;Mat dst;unsigned char * sImage, *dImage;int x, y;int nWidth, nHeight, step;orig = imread("test1.jpg", 1); /* 打开图像 */nWidth = orig.cols;nHeight = orig.rows;step = orig.cols * orig.channels();dst.create(Size(nWidth, nHeight), CV_8UC3);/* 创建目标图像 */sImage = (unsigned char*)malloc(sizeof(unsigned char)*(nHeight*nWidth * 3)); /* 创建2个图像buffer */dImage = (unsigned char*)malloc(sizeof(unsigned char)*(nHeight*nWidth * 3));/* 创建2个显示窗口，一个用于目标图像，一个用于源图像 */namedWindow("Original Image", WINDOW_AUTOSIZE);namedWindow("Result Image", WINDOW_AUTOSIZE);/* 取图像进行处理 */imshow("Original Image", orig);for (y = 0; y < orig.rows; y++){uchar* data = orig.ptr<uchar>(y);for (x = 0; x < orig.cols; x++){sImage[(y*nWidth + x) * 3 + 0] = data[3 * x + 0];sImage[(y*nWidth + x) * 3 + 1] = data[3 * x + 1];sImage[(y*nWidth + x) * 3 + 2] = data[3 * x + 2];}}/* 彩色图像增强 */MSRCR(sImage, nWidth, nHeight, orig.channels());for (y = 0; y < nHeight; y++){uchar* data = dst.ptr<uchar>(y);for (x = 0; x < nWidth; x++){data[3 * x + 0] = sImage[(y*nWidth + x) * 3 + 0];data[3 * x + 1] = sImage[(y*nWidth + x) * 3 + 1];data[3 * x + 2] = sImage[(y*nWidth + x) * 3 + 2];}}//显示处理图像imshow("Result Image", dst);waitKey(0);free(sImage);free(dImage);return 0;
}

下面的MSRCR是整个算法的和核心部分

/* 这个函数是算法的核心 */
/* (a) 在多个刻度处过滤，并将结果汇总 */
/* (b) 计算最终结果 */
void MSRCR(unsigned char * src, int width, int height, int bytes)
{int           scale, row, col;int           i, j;int           size;int           pos;int           channel;unsigned char *psrc = NULL;              /* SRC图的备份指针 */float         *dst = NULL;                     /* 算法的浮动缓冲区 */float         *pdst = NULL;                    /* 浮点缓冲区的备份指针 */float         *in, *out;int           channelsize;            /* 一个通道的浮动内存缓存 */float         weight;gauss3_coefs  coef;float         mean, var;float         mini, range, maxi;float         alpha;float         gain;float         offset;/* 分配算法所需的所有内存 */size = width * height * bytes;dst = (float *)malloc(size * sizeof(float));if (dst == NULL){printf("Failed to allocate memory");return;}memset(dst, 0, size * sizeof(float));channelsize = (width * height);in = (float *)malloc(channelsize * sizeof(float));if (in == NULL){free(dst);printf("Failed to allocate memory");return; /* 判断输入*/}out = (float *)malloc(channelsize * sizeof(float));if (out == NULL){free(in);free(dst);printf("Failed to allocate memory");return; /* 判断输出 */}/* 根据过滤器数量及其分布计算过滤尺度 */retinex_scales_distribution(RetinexScales, rvals.nscales, rvals.scales_mode, rvals.scale);/* 根据不同的尺度数进行过滤 *//* 根据一个特定的重量总结各种过滤器的结果（这里等同于所有过滤器）*/weight = 1.0f / (float)rvals.nscales;/* 根据所选的尺度，递归滤波算法需要不同的系数（~=高斯标准偏差）*/pos = 0;for (channel = 0; channel < 3; channel++){for (i = 0, pos = channel; i < channelsize; i++, pos += bytes){/* 0-255 => 1-256 */in[i] = (float)(src[pos] + 1.0);}//对原始图像进行每个尺度的高斯模糊for (scale = 0; scale < rvals.nscales; scale++){compute_coefs3(&coef, RetinexScales[scale]);/* 过滤（平滑）高斯递归 *//* 先筛选行 */for (row = 0; row < height; row++){pos = row * width;gausssmooth(in + pos, out + pos, width, 1, &coef);}memcpy(in, out, channelsize * sizeof(float));memset(out, 0, channelsize * sizeof(float));/* 过滤（平滑）高斯递归 *//* 选择列 */for (col = 0; col < width; col++){pos = col;gausssmooth(in + pos, out + pos, height, width, &coef);}/* 汇总过滤的值 *///对每个尺度下进行累加计算  Log[R(x,y)] =  Log[R(x,y)] + Weight(i)* ( Log[Ii(x,y)]-Log[Li(x,y)]);//其中Weight(i)表示每个尺度对应的权重，要求各尺度权重之和必须为1，经典的取值为等权重for (i = 0, pos = channel; i < channelsize; i++, pos += bytes){dst[pos] += weight * (float)(log(src[pos] + 1.0f) - log(out[i]));}}}free(in);free(out);/* 最终计算原始值和累积滤波器值 *//*   参数增益、alpha和offset是常量 *//* Ci(x,y)=log[a Ii(x,y)]-log[ Ei=1-s Ii(x,y)] */alpha = 128.0f;gain = 1.0f;offset = 0.0f;for (i = 0; i < size; i += bytes){float logl;psrc = src + i;pdst = dst + i;logl = (float)log((float)psrc[0] + (float)psrc[1] + (float)psrc[2] + 3.0f);pdst[0] = gain * ((float)(log(alpha * (psrc[0] + 1.0f)) - logl) * pdst[0]) + offset;pdst[1] = gain * ((float)(log(alpha * (psrc[1] + 1.0f)) - logl) * pdst[1]) + offset;pdst[2] = gain * ((float)(log(alpha * (psrc[2] + 1.0f)) - logl) * pdst[2]) + offset;}/* 根据一阶和二阶的统计数据调整颜色的动态 *//* 使用方差可以控制颜色的饱和度 */pdst = dst;//计算RGB各通道数据的均值Mean和均方差Varcompute_mean_var(pdst, &mean, &var, size, bytes);mini = mean - rvals.cvar*var;maxi = mean + rvals.cvar*var;range = maxi - mini;if (!range) range = 1.0;//对Log[R(x,y)]的每一个值，进行线性映射for (i = 0; i < size; i += bytes){psrc = src + i;pdst = dst + i;for (j = 0; j < 3; j++){float c = 255 * (pdst[j] - mini) / range;psrc[j] = (unsigned char)clip(c, 0, 255);}}free(dst);
}

计算期望分布的尺度

void retinex_scales_distribution(float* scales, int nscales, int mode, int s)
{if (nscales == 1){/* 一个参数选择中等尺度 */scales[0] = (float)s / 2;}else if (nscales == 2){/* 两个参数选择中和大尺度*/scales[0] = (float)s / 2;scales[1] = (float)s;}else{float size_step = (float)s / (float)nscales;int   i;/* 按照处理模式进行处理 */switch (mode){case RETINEX_UNIFORM:for (i = 0; i < nscales; ++i)scales[i] = 2.0f + (float)i * size_step;break;case RETINEX_LOW:size_step = (float)log(s - 2.0f) / (float)nscales;for (i = 0; i < nscales; ++i)scales[i] = 2.0f + (float)pow(10, (i * size_step) / log(10));break;case RETINEX_HIGH:size_step = (float)log(s - 2.0) / (float)nscales;for (i = 0; i < nscales; ++i)scales[i] = s - (float)pow(10, (i * size_step) / log(10));break;default:break;}}
}

计算一次平均值和方差

void compute_mean_var(float *src, float *mean, float *var, int size, int bytes)
{float vsquared;int i, j;float *psrc;vsquared = 0.0f;*mean = 0.0f;for (i = 0; i < size; i += bytes){psrc = src + i;for (j = 0; j < 3; j++){/* 平均值 */*mean += psrc[j];/* 方差 */vsquared += psrc[j] * psrc[j];}}*mean /= (float)size; /* mean */vsquared /= (float)size; /* mean (x^2) */*var = (vsquared - (*mean * *mean));*var = (float)sqrt(*var); /* var */
}

高斯模糊的快速计算

关于高斯模糊如果有不清楚的可以参考这篇文章

https://www.cnblogs.com/invisible2/p/9177018.html

void compute_coefs3(gauss3_coefs * c, float sigma)
{float q, q2, q3;q = 0;if (sigma >= 2.5f){q = 0.98711f * sigma - 0.96330f;}else if ((sigma >= 0.5f) && (sigma < 2.5f)){q = 3.97156f - 4.14554f * (float)sqrt((double)1 - 0.26891 * sigma);}else{q = 0.1147705018520355224609375f;}q2 = q * q;q3 = q * q2;c->b[0] = (1.57825f + (2.44413f*q) + (1.4281f *q2) + (0.422205f*q3));c->b[1] = ((2.44413f*q) + (2.85619f*q2) + (1.26661f *q3));c->b[2] = (-((1.4281f*q2) + (1.26661f *q3)));c->b[3] = ((0.422205f*q3));c->B = 1.0f - ((c->b[1] + c->b[2] + c->b[3]) / c->b[0]);c->sigma = sigma;c->N = 3;
}

高斯平滑

void gausssmooth(float *in, float *out, int size, int rowstride, gauss3_coefs *c)
{int i, n, bufsize;float *w1, *w2;/* forward pass */bufsize = size + 3;size -= 1;w1 = (float *)malloc(bufsize * sizeof(float));w2 = (float *)malloc(bufsize * sizeof(float));w1[0] = in[0];w1[1] = in[0];w1[2] = in[0];for (i = 0, n = 3; i <= size; i++, n++){w1[n] = (float)(c->B*in[i*rowstride] +((c->b[1] * w1[n - 1] +c->b[2] * w1[n - 2] +c->b[3] * w1[n - 3]) / c->b[0]));}/* backward pass */w2[size + 1] = w1[size + 3];w2[size + 2] = w1[size + 3];w2[size + 3] = w1[size + 3];for (i = size, n = i; i >= 0; i--, n--){w2[n] = out[i * rowstride] = (float)(c->B*w1[n] +((c->b[1] * w2[n + 1] +c->b[2] * w2[n + 2] +c->b[3] * w2[n + 3]) / c->b[0]));}free(w1);free(w2);
}

输入的原图

经过retinex算法处理的图，可以看到原图几乎什么都看不清，但是retinex算法可以将图片的颜色较好的还原出来。

代码已经传到了csdn上，下载链接是

https://download.csdn.net/download/weixin_42521239/11165325

没有积分的小伙伴可以在评论中留下你的邮箱，我发给你

图像增强--视网膜皮层Retinex算法（二）相关推荐

【图像增强】基于matlab双边滤波retinex算法暗光图像增强【含Matlab源码 2305期】
⛄一.简介 1 Retinex 1.1 理论 Retinex理论始于Land和McCann于20世纪60年代作出的一系列贡献,其基本思想是人感知到某点的颜色和亮度并不仅仅取决于该点进入人眼的绝对光线, ...
图像增强去雾之直方图均衡化/同态滤波/Retinex算法
图像增强去雾之直方图均衡化/同态滤波/Retinex算法最近撸了一发图像去雾的算法,主要举四个例子,分别用了全局直方图均衡化,局部直方图均衡化,同态滤波,Retinex增强算法.感兴趣的可以一起讨论 ...
Retinex算法解读
Retinex是一种常用的建立在科学实验和科学分析基础上的图像增强方法,它是Edwin.H.Land于1963年提出的.就跟Matlab是由Matrix和Laboratory合成的一样,Retinex ...
基于不均匀光照下的颜色校正——retinex算法，通态滤波算法
retinex算法原理及算法实现 Retinex是一种常用的建立在科学实验和科学分析基础上的图像增强方法,它是Edwin.H.Land于1963年提出的.就跟Matlab是由Matrix和Labora ...
Retinex、log对数变换、直方图均衡化区别，边缘增强Retinex算法与拉普拉斯算法联系、均衡化与亮度调节算法、大津阈值计算
1.其中Retinex算法具有的功能:动态范围压缩(即滤掉了低频部分,提取了高频).色调再现(即还有图像色彩):具有锐化.颜色恒常性.动态范围压缩大.色彩保真度高等特点. 从算法公式上的个人理 ...
使用pytorch从零开始实现YOLO-V3目标检测算法 (二)
原文:https://blog.csdn.net/u011520516/article/details/80212960 博客翻译这是从零开始实现YOLO v3检测器的教程的第2部分.在上一节中,我 ...
从零开始学数据结构和算法(二)线性表的链式存储结构
链表链式存储结构定义线性表的链式存储结构的特点是用一组任意的存储单元的存储线性表的数据元素,这组存储单元是可以连续的,也可以是不连续的. 种类结构图单链表应用:MessageQueue 插 ...
Unicode双向算法详解(bidi算法)(二)
作者:黄邦勇帅(原名:黄勇)2019-10-17 Unicode双向算法详解(bidi算法)(二) 本文为原创文章,转载请注明出处,或注明转载自"黄邦勇帅(原名:黄勇) 本文是对<C+ ...
票据ticket实现方式java代码_Java代码实践12306售票算法(二)
周五闲来无事,基于上一篇关于浅析12306售票算法(java版)理论,进行了java编码实践供各位读者参考(以下为相关代码的简单描述) 1.订票工具类 1.1初始化一列车厢的票据信息 /** * 生成 ...
matlab算法(二维傅立叶级数变换)
说明 Y = fft2(X) 使用快速傅里叶变换算法返回矩阵的二维傅里叶变换,这等同于计算 fft(fft(X).').'.如果 X 是一个多维数组,fft2 将采用高于 2 的每个维度的二维变换.输 ...

图像增强--视网膜皮层Retinex算法（二）

图像增强--视网膜皮层Retinex算法（二）相关推荐

最新文章

热门文章