自适应直方图均衡（CLAHE) 代码及详细注释【OpenCV】
【摘要】
效果图
标题
CLAHE简介
HE 直方图增强，大家都不陌生，是一种比较古老的对比度增强算法，它有两种变体：AHE 和 CLAHE；两者都是自适应的增强算法，功能差不多，但是前...
效果图
标题
CLAHE简介
HE 直方图增强，大家都不陌生，是一种比较古老的对比度增强算法，它有两种变体：AHE 和 CLAHE；两者都是自适应的增强算法，功能差不多，但是前者有一个很大的缺陷，就是有时候会过度放大图像中相同区域的噪声，为了解决这一问题，出现了 HE 的另一种改进算法，就是 CLAHE；CLAHE 是另外一种直方图均衡算法，CLAHE 和 AHE 的区别在于前者对区域对比度实行了限制，并且利用插值来加快计算。它能有效的增强或改善图像（局部）对比度，从而获取更多图像相关边缘信息有利于分割。还能够有效改善 AHE 中放大噪声的问题。另外，CLAHE 的有一个用途是被用来对图像去雾。
详细理论请参考博客
OpenCV源码的本地路径： %OPENCV%\opencv\sources\modules\imgproc\src\clahe.cpp
clahe.cpp
-
// ----------------------------------------------------------------------
-
// CLAHE
-
-
namespace
-
{
-
class CLAHE_CalcLut_Body : public cv::ParallelLoopBody
-
{
-
public:
-
CLAHE_CalcLut_Body(const cv::Mat& src, cv::Mat& lut, cv::Size tileSize, int tilesX, int clipLimit, float lutScale) :
-
src_(src), lut_(lut), tileSize_(tileSize), tilesX_(tilesX), clipLimit_(clipLimit), lutScale_(lutScale)
-
{
-
}
-
-
void operator ()(const cv::Range& range) const;
-
-
private:
-
cv::Mat src_;
-
mutable cv::Mat lut_;
-
-
cv::Size tileSize_;
-
int tilesX_;
-
int clipLimit_;
-
float lutScale_;
-
};
-
-
// 计算直方图查找表
-
void CLAHE_CalcLut_Body::operator ()(const cv::Range& range) const
-
{
-
const int histSize = 256;
-
-
uchar* tileLut = lut_.ptr(range.start);
-
const size_t lut_step = lut_.step; // size = tilesX_*tilesY_ * lut_step
-
-
// Range(0, tilesX_ * tilesY_)，全图被分为tilesX_*tiles_Y个块
-
for (int k = range.start; k < range.end; ++k, tileLut += lut_step)
-
{
-
// (tx, ty)表示当前所在是哪一块
-
// (0, 0) (1, 0)...(tilesX_-1, 0)
-
// (0, 1) (1, 1)...(tilesX_-1, 1)
-
// ...
-
// (0, tilesY_-1)... (tilesX_-1, tilesY_-1)
-
const int ty = k / tilesX_;
-
const int tx = k % tilesX_;
-
-
// retrieve tile submatrix
-
// 注意：tileSize.width表示分块的宽度，tileSize.height表示分块高度
-
cv::Rect tileROI;
-
tileROI.x = tx * tileSize_.width; // 换算为全局坐标
-
tileROI.y = ty * tileSize_.height;
-
tileROI.width = tileSize_.width;
-
tileROI.height = tileSize_.height;
-
-
const cv::Mat tile = src_(tileROI);
-
-
// calc histogram
-
int tileHist[histSize] = { 0, };
-
// 统计 ROI 的直方图
-
int height = tileROI.height;
-
const size_t sstep = tile.step;
-
for (const uchar* ptr = tile.ptr<uchar>(0); height--; ptr += sstep)
-
{
-
int x = 0;
-
for (; x <= tileROI.width - 4; x += 4)
-
{
-
int t0 = ptr[x], t1 = ptr[x + 1];
-
tileHist[t0]++; tileHist[t1]++;
-
t0 = ptr[x + 2]; t1 = ptr[x + 3];
-
tileHist[t0]++; tileHist[t1]++;
-
}
-
-
for (; x < tileROI.width; ++x)
-
tileHist[ptr[x]]++;
-
}
-
-
// clip histogram
-
if (clipLimit_ > 0)
-
{
-
// how many pixels were clipped
-
int clipped = 0;
-
for (int i = 0; i < histSize; ++i)
-
{
-
// 超过裁剪阈值
-
if (tileHist[i] > clipLimit_)
-
{
-
clipped += tileHist[i] - clipLimit_;
-
tileHist[i] = clipLimit_;
-
}
-
}
-
-
// redistribute clipped pixels
-
int redistBatch = clipped / histSize;
-
int residual = clipped - redistBatch * histSize;
-
-
// 平均分配裁剪的差值到所有直方图
-
for (int i = 0; i < histSize; ++i)
-
tileHist[i] += redistBatch;
-
-
// 处理差值的余数
-
for (int i = 0; i < residual; ++i)
-
tileHist[i]++;
-
}
-
-
// calc Lut
-
int sum = 0;
-
for (int i = 0; i < histSize; ++i)
-
{
-
// 累加直方图
-
sum += tileHist[i];
-
tileLut[i] = cv::saturate_cast<uchar>(sum * lutScale_); // static_cast<float>(histSize - 1) / tileSizeTotal
-
}
-
}
-
}
-
-
class CLAHE_Interpolation_Body : public cv::ParallelLoopBody
-
{
-
public:
-
CLAHE_Interpolation_Body(const cv::Mat& src, cv::Mat& dst, const cv::Mat& lut, cv::Size tileSize, int tilesX, int tilesY) :
-
src_(src), dst_(dst), lut_(lut), tileSize_(tileSize), tilesX_(tilesX), tilesY_(tilesY)
-
{
-
}
-
-
void operator ()(const cv::Range& range) const;
-
-
private:
-
cv::Mat src_;
-
mutable cv::Mat dst_;
-
cv::Mat lut_;
-
-
cv::Size tileSize_;
-
int tilesX_;
-
int tilesY_;
-
};
-
-
// 根据相邻4块的直方图插值
-
void CLAHE_Interpolation_Body::operator ()(const cv::Range& range) const
-
{
-
const size_t lut_step = lut_.step;
-
// Range(0, src.rows)
-
for (int y = range.start; y < range.end; ++y)
-
{
-
const uchar* srcRow = src_.ptr<uchar>(y);
-
uchar* dstRow = dst_.ptr<uchar>(y);
-
-
const float tyf = (static_cast<float>(y) / tileSize_.height) - 0.5f;
-
-
int ty1 = cvFloor(tyf);
-
int ty2 = ty1 + 1;
-
-
// 差值作为插值的比例
-
const float ya = tyf - ty1;
-
-
ty1 = std::max(ty1, 0);
-
ty2 = std::min(ty2, tilesY_ - 1);
-
-
const uchar* lutPlane1 = lut_.ptr(ty1 * tilesX_); // 当前块的直方图
-
const uchar* lutPlane2 = lut_.ptr(ty2 * tilesX_); // 向下一块的直方图
-
-
for (int x = 0; x < src_.cols; ++x)
-
{
-
const float txf = (static_cast<float>(x) / tileSize_.width) - 0.5f;
-
-
int tx1 = cvFloor(txf);
-
int tx2 = tx1 + 1;
-
-
// 差值作为插值的比例
-
const float xa = txf - tx1;
-
-
tx1 = std::max(tx1, 0);
-
tx2 = std::min(tx2, tilesX_ - 1);
-
-
// src_.ptr<uchar>(y)[x]
-
const int srcVal = srcRow[x];
-
-
// 索引 LUT
-
const size_t ind1 = tx1 * lut_step + srcVal;
-
const size_t ind2 = tx2 * lut_step + srcVal; // 向右一块的直方图
-
-
float res = 0;
-
// 根据直方图的值进行插值
-
// lut_.ptr(ty1 * tilesX_)[tx1 * lut_step + srcVa] => lut_[ty1][tx1][srcVal]
-
res += lutPlane1[ind1] * ((1.0f - xa) * (1.0f - ya));
-
res += lutPlane1[ind2] * ((xa) * (1.0f - ya));
-
res += lutPlane2[ind1] * ((1.0f - xa) * (ya));
-
res += lutPlane2[ind2] * ((xa) * (ya));
-
-
dstRow[x] = cv::saturate_cast<uchar>(res);
-
}
-
}
-
}
-
-
class CLAHE_Impl : public cv::CLAHE
-
{
-
public:
-
CLAHE_Impl(double clipLimit = 40.0, int tilesX = 8, int tilesY = 8);
-
-
cv::AlgorithmInfo* info() const; // Algorithm类工厂方法封装相关
-
-
void apply(cv::InputArray src, cv::OutputArray dst);
-
-
void setClipLimit(double clipLimit);
-
double getClipLimit() const;
-
-
void setTilesGridSize(cv::Size tileGridSize);
-
cv::Size getTilesGridSize() const;
-
-
void collectGarbage();
-
-
private:
-
double clipLimit_;
-
int tilesX_;
-
int tilesY_;
-
-
cv::Mat srcExt_;
-
cv::Mat lut_;
-
};
-
-
CLAHE_Impl::CLAHE_Impl(double clipLimit, int tilesX, int tilesY) :
-
clipLimit_(clipLimit), tilesX_(tilesX), tilesY_(tilesY)
-
{
-
}
-
// Algorithm类工厂方法封装相关
-
//CV_INIT_ALGORITHM(CLAHE_Impl, "CLAHE",
-
// obj.info()->addParam(obj, "clipLimit", obj.clipLimit_);
-
//obj.info()->addParam(obj, "tilesX", obj.tilesX_);
-
//obj.info()->addParam(obj, "tilesY", obj.tilesY_))
-
-
void CLAHE_Impl::apply(cv::InputArray _src, cv::OutputArray _dst)
-
{
-
cv::Mat src = _src.getMat();
-
-
CV_Assert(src.type() == CV_8UC1);
-
-
_dst.create(src.size(), src.type());
-
cv::Mat dst = _dst.getMat();
-
-
const int histSize = 256;
-
// 准备 LUT，tilesX_*tilesY_个块，每个块都有256个柱子的直方图
-
lut_.create(tilesX_ * tilesY_, histSize, CV_8UC1);
-
-
cv::Size tileSize;
-
cv::Mat srcForLut;
-
-
// 如果分块刚好（整除）
-
if (src.cols % tilesX_ == 0 && src.rows % tilesY_ == 0)
-
{
-
tileSize = cv::Size(src.cols / tilesX_, src.rows / tilesY_);
-
srcForLut = src;
-
}
-
// 否则对原图进行扩充
-
else
-
{
-
cv::copyMakeBorder(src, srcExt_, 0, tilesY_ - (src.rows % tilesY_), 0, tilesX_ - (src.cols % tilesX_), cv::BORDER_REFLECT_101);
-
-
tileSize = cv::Size(srcExt_.cols / tilesX_, srcExt_.rows / tilesY_);
-
srcForLut = srcExt_;
-
}
-
-
const int tileSizeTotal = tileSize.area();
-
const float lutScale = static_cast<float>(histSize - 1) / tileSizeTotal; // △
-
-
// 计算实际的clipLimit
-
int clipLimit = 0;
-
if (clipLimit_ > 0.0)
-
{
-
clipLimit = static_cast<int>(clipLimit_ * tileSizeTotal / histSize);
-
clipLimit = std::max(clipLimit, 1);
-
}
-
-
// 分块并行计算： LUT
-
CLAHE_CalcLut_Body calcLutBody(srcForLut, lut_, tileSize, tilesX_, clipLimit, lutScale);
-
cv::parallel_for_(cv::Range(0, tilesX_ * tilesY_), calcLutBody);
-
-
// 分块并行计算： 根据直方图插值
-
CLAHE_Interpolation_Body interpolationBody(src, dst, lut_, tileSize, tilesX_, tilesY_);
-
cv::parallel_for_(cv::Range(0, src.rows), interpolationBody);
-
}
-
-
void CLAHE_Impl::setClipLimit(double clipLimit)
-
{
-
clipLimit_ = clipLimit;
-
}
-
-
double CLAHE_Impl::getClipLimit() const
-
{
-
return clipLimit_;
-
}
-
-
void CLAHE_Impl::setTilesGridSize(cv::Size tileGridSize)
-
{
-
tilesX_ = tileGridSize.width;
-
tilesY_ = tileGridSize.height;
-
}
-
-
cv::Size CLAHE_Impl::getTilesGridSize() const
-
{
-
return cv::Size(tilesX_, tilesY_);
-
}
-
-
void CLAHE_Impl::collectGarbage()
-
{
-
srcExt_.release();
-
lut_.release();
-
}
-
}
-
-
cv::Ptr<cv::CLAHE> cv::createCLAHE(double clipLimit, cv::Size tileGridSize)
-
{
-
return new CLAHE_Impl(clipLimit, tileGridSize.width, tileGridSize.height);
-
}
main.cpp
-
int main(int argc, char** argv)
-
{
-
cv::Mat inp_img = cv::imread("D:/Pictures/beard.jpg");
-
if (!inp_img.data) {
-
cout << "Something Wrong";
-
return -1;
-
}
-
namedWindow("Input Image", CV_WINDOW_AUTOSIZE);
-
cv::imshow("Input Image", inp_img);
-
-
cv::Mat clahe_img;
-
cv::cvtColor(inp_img, clahe_img, CV_BGR2Lab);
-
std::vector<cv::Mat> channels(3);
-
cv::split(clahe_img, channels);
-
-
cv::Ptr<cv::CLAHE> clahe = cv::createCLAHE();
-
// 直方图的柱子高度大于计算后的ClipLimit的部分被裁剪掉，然后将其平均分配给整张直方图
-
// 从而提升整个图像
-
clahe->setClipLimit(4.); // (int)(4.*(8*8)/256)
-
//clahe->setTilesGridSize(Size(8, 8)); // 将图像分为8*8块
-
cv::Mat dst;
-
clahe->apply(channels[0], dst);
-
dst.copyTo(channels[0]);
-
cv::merge(channels, clahe_img);
-
-
cv::Mat image_clahe;
-
cv::cvtColor(clahe_img, image_clahe, CV_Lab2BGR);
-
-
//cout << cvFloor(-1.5) << endl;
-
-
namedWindow("CLAHE Image", CV_WINDOW_AUTOSIZE);
-
cv::imshow("CLAHE Image", image_clahe);
-
imwrite("out.jpg", image_clahe);
-
cv::waitKey(0);
-
destroyAllWindows();
-
-
return 0;
-
}
注意：cv::ParallelLoopBody 位于 %OpenCV%\opencv\sources\modules\core\src\parallel.cpp
延伸阅读：Algorithm类工厂方法封装相关
文章来源: panda1234lee.blog.csdn.net，作者：panda1234lee，版权归原作者所有，如需转载，请联系作者。
原文链接：panda1234lee.blog.csdn.net/article/details/52852765
【版权声明】本文为华为云社区用户转载文章，如果您发现本社区中有涉嫌抄袭的内容，欢迎发送邮件进行举报，并提供相关证据，一经查实，本社区将立刻删除涉嫌侵权内容，举报邮箱：
cloudbbs@huaweicloud.com
- 点赞
- 收藏
- 关注作者
评论（0）