圖像處理——（源）中值濾波（MediFilter）函數編程實現

 

https://blog.csdn.net/tengfei461807914/article/details/83626123

中值濾波是一種非線性濾波，在處理脈沖噪聲以及椒鹽噪聲時效果極佳，能夠有效的保護好圖像的邊緣信息。

中值濾波的處理思路很簡單，取卷積核當中所覆蓋像素中的中值作為錨點的像素值即可。

如果按照遍歷所有像素，再對卷積核中的像素排序取中值，那么時間復雜度會很高，需要對中值濾波進行改進。

中值濾波的改進實際上很是很好想的，無非就是一個滑動窗口取中值的問題，每次向右滑動的過程中等於在窗口中新添加添加一列窗口像素，同時減去一列窗口像素，考慮維護這個窗口中的像素信息變化即可。

這里面使用huang算法，該思路是這個人提出來的，算法思路如下：

在計算中值的辦法中，不使用排序，而是使用像素直方圖，也就是記錄像素值的哈希。首先設定閾值threshold，這個threshold就是窗口的中心位置，即ksize×ksize/2+1，kisze為窗口尺寸。

每次在計算中值的過程中，從小到大累加像素直方圖的值，如果該值大於等於，此時對應的像素值就是中值了。

例如ksize=3的窗口如下：
⎡⎣⎢122235154⎤⎦⎥(3) \left[\begin{matrix}1 & 2 & 1 \\2 & 3 & 5 \\2 & 5 & 4\end{matrix}\right] \tag{3}
⎣
⎡
​

1
2
2
​

2
3
5
​

1
5
4
​

⎦
⎤
​
(3)

對該窗口中的值計算像素直方圖如下，threshold=3×3/2+1=5
1:2（表示像素值為1的有2個）
2:3
3:1
4:1
5:2

因為2+3≥5，所以中值為2

每次滑動窗口的過程中，如果窗口在第一列，那么直接計算直方圖。否則向右移動，在直方圖中減去左側離開窗口中像素，在右側添加進入窗口中的像素。

此外，也可以讓窗口按照蛇形來移動，這樣也會避免每次窗口在第一列時需要重新計算的問題。

void AddSultPapperNoise(const Mat &src, Mat &dst,int num)//添加椒鹽噪聲
{
    dst = src.clone();
    uchar *pd=dst.data;
    int row, col, cha;
    srand((unsigned)time(NULL));
    while (num--)
    {
        row = rand() % dst.rows;
        col = rand() % dst.cols;
        cha = rand() % dst.channels();
        pd[(row*dst.cols + col)*dst.channels() + cha] = 0;
    }
}


int GetMediValue(const int histogram[], int thresh)//計算中值
{
    int sum = 0;
    for (int i = 0; i < (1 << 16); i++)
    {
        sum += histogram[i];
        if (sum >= thresh)
            return i;
    }
    return (1 << 16);
}

void MyFastMediFilter(const Mat &src, Mat &dst, int ksize)
{
    CV_Assert(ksize % 2 == 1);

    Mat tmp;
    int len = ksize / 2;
    tmp.create(Size(src.cols + len, src.rows + len), src.type());//添加邊框
    dst.create(Size(src.cols, src.rows), src.type());


    int channel = src.channels();
    uchar *ps = src.data;
    uchar *pt = tmp.data;
    for (int row = 0; row < tmp.rows; row++)//添加邊框的過程
    {
        for (int col = 0; col < tmp.cols; col++)
        {
            for (int c = 0; c < channel; c++)
            {
                if (row >= len && row < tmp.rows - len && col >= len && col < tmp.cols - len)
                    pt[(tmp.cols * row + col) * channel + c] = ps[(src.cols * (row - len) + col - len) * channel + c];
                else
                    pt[(tmp.cols * row + col) * channel + c] = 0;
            }
        }
    }
    int Hist[(1 << 16)] = { 0 };
    uchar *pd = dst.data;
    ushort val = 0;
    pt = tmp.data;
    for (int c = 0; c < channel; c++)//每個通道單獨計算
    {
        for (int row = len; row < tmp.rows - len; row++)
        {
            for (int col = len; col < tmp.cols - len; col++)
            {
                
                if (col == len)
                {
                    memset(Hist, 0, sizeof(Hist));
                    for (int x = -len; x <= len; x++)
                    {
                        for (int y = -len; y <= len; y++)
                        {
                            val = pt[((row + x) * tmp.cols + col + y) * channel + c];
                            Hist[val]++;
                        }
                    }
                }
                else
                {
                    int L = col - len - 1;
                    int R = col + len;
                    for (int y = -len; y <= len; y++)
                    {
                        int leftInd = ((row + y) * tmp.cols + L) * channel + c;
                        int rightInd = ((row + y) * tmp.cols + R) * channel + c;
                        Hist[pt[leftInd]]--;
                        Hist[pt[rightInd]]++;
                    }
                }
                val = GetMediValue(Hist, ksize*ksize / 2 + 1);
                pd[(dst.cols * (row - len) + col - len) * channel + c] = val;
                
            }
        }
    }
}