自适应直方图均衡(CLAHE) 代码及详细注释【OpenCV】

本文涉及的产品
云数据库 Tair(兼容Redis),内存型 2GB
Redis 开源版,标准版 2GB
推荐场景:
搭建游戏排行榜
简介: 效果图

效果图

image.pngimage.png

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

相关实践学习
基于Redis实现在线游戏积分排行榜
本场景将介绍如何基于Redis数据库实现在线游戏中的游戏玩家积分排行榜功能。
云数据库 Redis 版使用教程
云数据库Redis版是兼容Redis协议标准的、提供持久化的内存数据库服务,基于高可靠双机热备架构及可无缝扩展的集群架构,满足高读写性能场景及容量需弹性变配的业务需求。 产品详情:https://www.aliyun.com/product/kvstore &nbsp; &nbsp; ------------------------------------------------------------------------- 阿里云数据库体验:数据库上云实战 开发者云会免费提供一台带自建MySQL的源数据库&nbsp;ECS 实例和一台目标数据库&nbsp;RDS实例。跟着指引,您可以一步步实现将ECS自建数据库迁移到目标数据库RDS。 点击下方链接,领取免费ECS&amp;RDS资源,30分钟完成数据库上云实战!https://developer.aliyun.com/adc/scenario/51eefbd1894e42f6bb9acacadd3f9121?spm=a2c6h.13788135.J_3257954370.9.4ba85f24utseFl
目录
相关文章
|
6月前
|
传感器 C++ 计算机视觉
【opencv3】详述PnP测距完整流程(附C++代码)
【opencv3】详述PnP测距完整流程(附C++代码)
|
6月前
|
机器学习/深度学习 算法 数据可视化
计算机视觉+深度学习+机器学习+opencv+目标检测跟踪+一站式学习(代码+视频+PPT)-2
计算机视觉+深度学习+机器学习+opencv+目标检测跟踪+一站式学习(代码+视频+PPT)
|
15天前
|
机器学习/深度学习 监控 算法
基于计算机视觉(opencv)的运动计数(运动辅助)系统-源码+注释+报告
基于计算机视觉(opencv)的运动计数(运动辅助)系统-源码+注释+报告
30 3
|
5月前
|
算法 开发工具 计算机视觉
【零代码研发】OpenCV实验大师工作流引擎C++ SDK演示
【零代码研发】OpenCV实验大师工作流引擎C++ SDK演示
75 1
|
3月前
|
算法 计算机视觉 Python
python利用opencv进行相机标定获取参数,并根据畸变参数修正图像附有全部代码(流畅无痛版)
该文章详细介绍了使用Python和OpenCV进行相机标定以获取畸变参数,并提供了修正图像畸变的全部代码,包括生成棋盘图、拍摄标定图像、标定过程和畸变矫正等步骤。
python利用opencv进行相机标定获取参数,并根据畸变参数修正图像附有全部代码(流畅无痛版)
|
3月前
|
计算机视觉 Python
opencv在pycharm不能自动补全代码
opencv在pycharm不能自动补全代码
30 0
|
6月前
|
机器学习/深度学习 Ubuntu Linux
计算机视觉+深度学习+机器学习+opencv+目标检测跟踪+一站式学习(代码+视频+PPT)-1
计算机视觉+深度学习+机器学习+opencv+目标检测跟踪+一站式学习(代码+视频+PPT)
|
5月前
|
监控 安全 计算机视觉
实战 | 18行代码轻松实现人脸实时检测【附完整代码与源码详解】Opencv、人脸检测
实战 | 18行代码轻松实现人脸实时检测【附完整代码与源码详解】Opencv、人脸检测
|
5月前
|
并行计算 IDE 开发工具
【竹篮打水】OpenCV4.x 中新增并行代码执行演示
【竹篮打水】OpenCV4.x 中新增并行代码执行演示
45 0
|
6月前
|
算法 API 计算机视觉
基于opencv的大米计数统计(详细处理流程+代码)
基于opencv的大米计数统计(详细处理流程+代码)
基于opencv的大米计数统计(详细处理流程+代码)