前言

之前做的一些滤波，例如高斯，双边，均值，导向滤波均是在时域下做的滤波。而同态滤波是在频域下来做滤波，用于改善动态范围很大但是暗区细节又不清楚的图像。

原理

对于一般的图像，有这样一个模型 $f (x, y) = i (x, y) r (x, y)$ 。其中f(x, y)代表原图像，i(x,y)代表照射强度，r(x,y)代表反射强度。一般的图像光照是均匀变化的，所以i应该是低频分量，而不同的物体对光的反射是具有突变型的，所以r(x,y)是高频分量。现在同时对上式两边同时取对数，在做Fourier变化，得到线性组合的频率域： $l n f (x, y) = l n i (x, y) + l n r (x, y)$ $F F T (l n f (x, y)) = F F T (l n i (x, y)) + F F T (l n r (x, y))$ ，我们希望对低频能量进行压制，这样就降低了动态范围，而要对高频分量进行提高，这样会增强图像的对比度。因此，同态滤波器的传递函数一般在低频部分小于1，高频部分大于1。

算法过程

对原图像f(x,y)取对数，目的是使得图像模型中的乘法运算转化为简单的加法运算： $z (x, y) = l n f (x, y) = l n i (x, y) + l n r (x, y)$
再对对数函数做傅里叶变换，将图像转换到频域
$F (z (x, y)) = F [l n i (x, y)] + F [l n r (x, y)]$ 即 $Z = I + R$
选择适当的传递函数，压缩照射分量 $i (x, y)$ 的变化范围，削弱 $I (u, v)$ ，增强反射分量 $r (x, y)$ 的对比度，所以总结起来就是这个滤波器需要对低频能量进行压制，以降低动态范围，同时要对高频进行提高，以增强图像对比度。所以同态滤波器的传递函数如下：
假设用一个同态滤波器函数 $H (u, v)$ 来处理原图像f(x,y)的对数的傅里叶变换 $Z (u, v)$ 得到： $S (u, v) = H (u, v) Z (u, v) = H (u, v) I (u, v) + H (u, v) R (u, v)$ ，逆变换到时域得： $s (x, y) = F^{- 1} (S (u, v))$ ，所以 $f^{'} (x, y) = e x p (s (x, y))$ ， $f^{'} (x, y)$ 代表滤波后的图像。整个算法的过程可以用下面的流程图来表示：

C++ opencv实现

传统的同态滤波器有高斯同态滤波 $H (u, v)_{1} = (R h - R l) [1 - e x p (- c (D (u, v) / D_{0})^{2 n})] + R l$ ，巴特沃斯同态滤波器 $H (u, v)_{2} = (R h - R l) [1 / (1 + [D_{0} / c D (u, v)]^{2 n})] + R I$ ,指数型同态滤波器： $H (u, v)_{3} = (R h - R l) e x p (- c [D_{0} / D (u, v)])^{n} + R I$ 。其中 $D (u, v) = [(u - u_{0})^{2} + (v - v_{0})^{2}]^{\frac{1}{2}}$ ,Rh,RI分别为高频和低频增益，D0为截止频率，c是控制斜面锐化的常数。这里实现了第一种滤波器，也就是高斯滤波。

Mat HomoFilter(cv::Mat src){
    src.convertTo(src, CV_64FC1);
    int rows = src.rows;
    int cols = src.cols;
    int m = rows % 2 == 1 ? rows + 1 : rows;
    int n = cols % 2 == 1 ? cols + 1 : cols;
    copyMakeBorder(src, src, 0, m - rows, 0, n - cols, BORDER_CONSTANT, Scalar::all(0));
    rows = src.rows;
    cols = src.cols;
    Mat dst(rows, cols, CV_64FC1);
    //1. ln
    for(int i = 0; i < rows; i++){
        double *srcdata = src.ptr<double>(i);
        double *logdata = src.ptr<double>(i);
        for(int j = 0; j < cols; j++){
            logdata[j] = log(srcdata[j] + 0.0001);
        }
    }
    //2. dct
    Mat mat_dct = Mat::zeros(rows, cols, CV_64FC1);
    dct(src, mat_dct);
    //3. 高斯同态滤波器
    Mat H_u_v;
    double gammaH = 1.5;
    double gammaL = 0.5;
    double C = 1;
    double  d0 = (src.rows / 2) * (src.rows / 2) + (src.cols / 2) * (src.cols / 2);
    double  d2 = 0;
    H_u_v = Mat::zeros(rows, cols, CV_64FC1);
    for(int i = 0; i < rows; i++){
        double * dataH_u_v = H_u_v.ptr<double>(i);
        for(int j = 0; j < cols; j++){
            d2 = pow(i, 2.0) + pow(j, 2.0);
            dataH_u_v[j] = (gammaH - gammaL) * (1 - exp(-C * d2 / d0)) + gammaL;
        }
    }
    H_u_v.ptr<double>(0)[0] = 1.1;
    mat_dct = mat_dct.mul(H_u_v);
    //4. idct
    idct(mat_dct, dst);
    //exp
    for(int i = 0; i < rows; i++){
        double  *srcdata = dst.ptr<double>(i);
        double *dstdata = dst.ptr<double>(i);
        for(int j = 0; j < cols; j++){
            dstdata[j] = exp(srcdata[j]);
        }
    }
    dst.convertTo(dst, CV_8UC1);
    return dst;
}

主函数：

int main(){
    //
    Mat src = cv::imread("./1.png");
   imshow("origin", src);
    int originrows = src.rows;
    int origincols = src.cols;
    Mat dst(src.rows, src.cols, CV_8UC3);
    cvtColor(src, src, COLOR_BGR2YUV);
    vector <Mat> yuv;
    split(src, yuv);
    Mat nowY = yuv[0];
    Mat newY = HomoFilter(nowY);
    Mat tempY(originrows, origincols, CV_8UC1);
    for(int i = 0; i < originrows; i++){
        for(int j = 0; j < origincols; j++){
            tempY.at<uchar>(i, j) = newY.at<uchar>(i, j);
        }
    }
    yuv[0] = tempY;
    merge(yuv, dst);
    cvtColor(dst, dst, COLOR_YUV2BGR);
    imshow("result", dst);
    waitKey(0);
}

效果：

可以看到图片的整体亮度被拉到人眼可以接受的范围里面了，但是缺点是有些地方出现了白点。再看一组效果：

我的这个代码效果很不好，如果真的要使用同态滤波需要认真取做研究，可以参考https://github.com/lilingyu/homofilter

同态滤波 原理及C++实现

前言

原理

算法过程

C++ opencv实现

同态滤波原理及C++实现