效果图
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类工厂方法封装相关
