ACE(自动色彩均衡算法)C++代码实现

参考内容:

python实现:

https://blog.csdn.net/zmshy2128/article/details/53470357?utm_medium=distribute.pc_relevant.none-task-blog-baidulandingword-2&spm=1001.2101.3001.4242

C++实现:

https://blog.csdn.net/just_sort/article/details/85237711

 

(该版本写的很清楚,但速度较慢,所以个人改写了部分,时间花费大约缩短三分之二)

主要起到提速作用的是采用Mat的指针寻址.

#include <stdio.h>
#include <iostream>
#include <immintrin.h>
#include <opencv2/opencv.hpp>
#include <opencv2/core/core.hpp>
#include <opencv2/ml/ml.hpp>
#include "opencv2/highgui/highgui.hpp"
#include <time.h>
using namespace cv;
using namespace cv::ml;
using namespace std;


namespace ACE {
    Mat stretchImage(Mat src) {
        int row = src.rows;
        int col = src.cols;
        Mat dst(row, col, CV_64FC1);
        double MaxValue = *max_element(src.begin<double>(), src.end<double>());
        double MinValue = *min_element(src.begin<double>(), src.end<double>());
        for (int i = 0; i < row; i++) {
            double *ptr_dst=dst.ptr<double>(i);
            double *ptr_src=src.ptr<double>(i);
            for (int j = 0; j < col; j++) {
                double src_data = *ptr_src++;
                *ptr_dst++=(src_data- MinValue) / (MaxValue - MinValue);
            }
        }
        return dst;
    }

    Mat getPara(int radius) {
        int size = radius * 2 + 1;
        Mat dst(size, size, CV_64FC1);
        for (int i = -radius; i <= radius; i++) {
            for (int j = -radius; j <= radius; j++) {
                if (i == 0 && j == 0) {
                    dst.at<double>(i + radius, j + radius) = 0;
                }
                else {
                    dst.at<double>(i + radius, j + radius) = 1.0 / sqrt(i * i + j * j);
                }
            }
        }
        double sum = 0;
        for (int i = 0; i < size; i++) {
            for (int j = 0; j < size; j++) {
                sum += dst.at<double>(i, j);
            }
        }
        for (int i = 0; i < size; i++) {
            for (int j = 0; j < size; j++) {
                dst.at<double>(i, j) = dst.at<double>(i, j) / sum;
            }
        }
        return dst;
    }

    Mat NormalACE(Mat src, int ratio, int radius) {
        Mat para = getPara(radius);
        int row = src.rows;
        int col = src.cols;
        int size = 2 * radius + 1;
        Mat Z(row + 2 * radius, col + 2 * radius, CV_64FC1);
        for (int i = 0; i < Z.rows; i++) {
            double *ptr_Z=Z.ptr<double>(i);
            double *ptr_src=src.ptr<double>(i-radius);
            for (int j = 0; j < Z.cols; j++) {
                if((i - radius >= 0) && (i - radius < row) && (j - radius >= 0) && (j - radius < col)) {
                    *ptr_Z++=ptr_src[j-radius];
                    //Z.at<double>(i, j) = src.at<double>(i - radius, j - radius);
                }
                else {
                    *ptr_Z++=0;
                    //Z.at<double>(i, j) = 0;
                }
            }
        }
        Mat dst=Mat::zeros(row,col,CV_64FC1);
        for (int i = 0; i < size; i++) {
            double *p_para=para.ptr<double>(i);
            for (int j = 0; j < size; j++) {
                double para_data=*p_para;
                if (*p_para++ == 0) continue;
                for (int x = 0; x < row; x++) {
                    double *p_src=src.ptr<double>(x);
                    double *p_dst=dst.ptr<double>(x);
                    double *p_Z=Z.ptr<double>(x+i,j);
                    for (int y = 0; y < col; y++) {
                        double scr_at_xy = *p_src++;
                        double sub = scr_at_xy -*p_Z++;// Z.at<double>(x + i, y + j);
                        //double sub = src.at<double>(x, y) - Z.at<double>(x + i, y + j);
                        double tmp = sub * ratio;
                        if (tmp > 1.0) tmp = 1.0;
                        if (tmp < -1.0) tmp = -1.0;
                        *(p_dst++)+=tmp * para_data;
                    }
                }
            }
        }

        return dst;
    }

    Mat FastACE(Mat src, int ratio, int radius) {
        int row = src.rows;
        int col = src.cols;
        if (min(row, col) <= 2) {
            Mat dst(row, col, CV_64FC1);
            for (int i = 0; i < row; i++) {
                for (int j = 0; j < col; j++) {
                    dst.at<double>(i, j) = 0.5;
                }
            }
            return dst;
        }

        Mat Rs((row + 1) / 2, (col + 1) / 2, CV_64FC1);
        resize(src, Rs, Size((col + 1) / 2, (row + 1) / 2));
        Mat Rf= FastACE(Rs, ratio, radius);
        resize(Rf, Rf, Size(col, row));
        resize(Rs, Rs, Size(col, row));
        Mat dst(row, col, CV_64FC1);
        Mat dst1 = NormalACE(src, ratio, radius);
        Mat dst2 = NormalACE(Rs, ratio, radius);

        dst=Rf+dst1-dst2;
        return dst;
    }

    Mat getACE(Mat src, int ratio, int radius) {
        int row = src.rows;
        int col = src.cols;
        vector <Mat> v;
        split(src, v);
        v[0].convertTo(v[0], CV_64FC1);
        v[1].convertTo(v[1], CV_64FC1);
        v[2].convertTo(v[2], CV_64FC1);
        Mat src1(row, col, CV_64FC1);
        Mat src2(row, col, CV_64FC1);
        Mat src3(row, col, CV_64FC1);

        for (int i = 0; i < row; i++) {
            for (int j = 0; j < col; j++) {
                src1.at<double>(i, j) = 1.0 * src.at<Vec3b>(i, j)[0] / 255.0;
                src2.at<double>(i, j) = 1.0 * src.at<Vec3b>(i, j)[1] / 255.0;
                src3.at<double>(i, j) = 1.0 * src.at<Vec3b>(i, j)[2] / 255.0;
            }
        }
        src1 = stretchImage(FastACE(src1, ratio, radius));
        src2 = stretchImage(FastACE(src2, ratio, radius));
        src3 = stretchImage(FastACE(src3, ratio, radius));
        Mat dst1(row, col, CV_8UC1);
        Mat dst2(row, col, CV_8UC1);
        Mat dst3(row, col, CV_8UC1);
        for (int i = 0; i < row; i++) {
            for (int j = 0; j < col; j++) {
                dst1.at<uchar>(i, j) = (int)(src1.at<double>(i, j) * 255);
                if (dst1.at<uchar>(i, j) > 255) dst1.at<uchar>(i, j) = 255;
                else if (dst1.at<uchar>(i, j) < 0) dst1.at<uchar>(i, j) = 0;
                dst2.at<uchar>(i, j) = (int)(src2.at<double>(i, j) * 255);
                if (dst2.at<uchar>(i, j) > 255) dst2.at<uchar>(i, j) = 255;
                else if (dst2.at<uchar>(i, j) < 0) dst2.at<uchar>(i, j) = 0;
                dst3.at<uchar>(i, j) = (int)(src3.at<double>(i, j) * 255);
                if (dst3.at<uchar>(i, j) > 255) dst3.at<uchar>(i, j) = 255;
                else if (dst3.at<uchar>(i, j) < 0) dst3.at<uchar>(i, j) = 0;
            }
        }

        vector <Mat> out;
        out.push_back(dst1);
        out.push_back(dst2);
        out.push_back(dst3);
        Mat dst;
        merge(out, dst);
        return dst;
    }
}

using namespace ACE;

int main() {
    Mat src = imread("./test2.png");
    while(1)
    {
        clock_t start_t = clock();  //时间起始
        Mat dst = getACE(src, 15, 1);
        clock_t end_t   = clock(); //时间测试结束
        std::cout<<"耗时:"<<double(end_t - start_t) / CLOCKS_PER_SEC <<std::endl;
        imshow("origin", src);
        imshow("result", dst);
        int ck=waitKey(0);
        if (ck==27)
            break;
    }
    cv::destroyAllWindows();
}