直方圖實現快速中值濾波

本文轉載自查看原文 2017-12-05 22:42 5120 OpenCV/ 圖像算法/ 中值濾波/ opencv

中值濾波能夠有效去除圖像中的異常點，具有去除圖像噪聲的作用。傳統中值濾波的算法一般都是在圖像中建立窗口，然后對窗口內的所有像素值進行排序，選擇排序后的中間值作為窗口中心像素濾波后的值。由於這個做法在每個像素點處都要建立窗口並排序，非常耗時，尤其是有大量的冗余計算。如下圖：

黃色區域+中間粉色區域是第一個像素為中心建立的濾波窗口，粉色區域+右邊藍色區域為同一行第二個像素為中心建立的濾波窗口。傳統做法對每一個窗口內所有像素排序，那么中間粉色區域的像素就會被重復的做排序運算，存在大量冗余計算。如果窗口移動一個像素的時候只用考慮去除左邊一列內（黃色區域）的像素，並且加上右邊一列（藍色區域）的像素，運算會大大簡化。這樣的操作可以使用直方圖來實現。

一、直方圖實現快速中值濾波算法流程：

1.讀取圖像I，並且設定濾波窗口大小(winX*winY)，一般winX=winY，奇數。

2.設定中值濾波直方圖中的閾值，Thresh=(winX*winY)/2 +1;

3.如果要考慮邊界情況，可以先對原圖像進行擴展，左、右邊界分別擴展winX/2個像素，上下邊界分別擴展winY/2個像素。

4.逐行遍歷圖像像素，以第一行為例：先取第一行第一個要處理的像素（窗口中心像素），建立濾波窗口，提取窗口內所有像素值（N=winX*winY個），獲取N個像素的直方圖Hist。從左到右累加直方圖中每個灰度層級下像素點個數，記為sumCnt，直到sumCnt>=Thresh，這時的灰度值就是當前窗口內所有像素值的中值MediaValue。將MediaValue值賦值給窗口中心像素，表明第一個像素中值濾波完成。

5.此時窗口需要向右移動一個像素，開始濾波第二個像素，並且更新直方圖。以第二個像素為窗口中心建立濾波窗口，從前一個窗口的灰度直方圖Hist中減去窗口中最左側的一列像素值的灰度個數，然后加上窗口最右側一列像素值的灰度個數。完成直方圖的更新。

6.直方圖更新后，sumCnt值有三種變化可能：（1）減小（2）維持不變（3）增大。這三種情況與減去與加入的像素值灰度有關。此時為了求得新的中值，需要不斷調整sumCnt與Thresh之間的關系。

（1）如果sumCnt值小於Thresh：說明中值在直方圖當前灰度層級的右邊，sumCnt就依次向右加上一個灰度層級中灰度值個數，直到滿足sumCnt>=Thresh為止。記錄此時的灰度層級代表的灰度值，更新MediaValue，作為第二個像素的濾波后的值。

（2）維持不變：說明MediaValue值不變，直接作為第二個像素濾波后的值。

（3）如果sumCnt值大於Thresh：說明中值在直方圖當前灰度層級的左邊，sumCnt就依次向左減去一個灰度層級中灰度值個數，直到滿足sumCnt<=Thresh為止。記錄此時的灰度層級代表的灰度值，更新MediaValue值，作為第二個像素的濾波后的值。

7.窗口逐行依次滑動，求得整幅圖像的中值濾波結果。

二、濾波結果

以下圖手機拍攝的moon.jpg為例：

OpenCV中值濾波結果：

直方圖快速濾波結果：

完整代碼（兩種實現，原理一樣）如下：（博主偷懶沒有提前做邊界擴展，而是直接保留了四個邊界的像素值，邊界擴展也很容易實現，不再贅述）

Code01:

  1 #include <opencv2\opencv.hpp>
  2 #include <iostream>
  3 #include <string>
  4 
  5 using namespace cv;
  6 using namespace std;
  7 
  8 //計算亮度中值和灰度<=中值的像素點個數
  9 int calMediaValue(const int histogram[], int thresh)
 10 {
 11     int sum = 0;
 12     for (int i = 0; i < 256; ++i)
 13     {
 14         sum += histogram[i];
 15         if (sum>= thresh)
 16         {
 17             return i;
 18         }
 19     }
 20     return 255;
 21 }
 22 
 23 void fastMedianBlur(const Mat &srcImg, Mat &dstImg, int diameter)
 24 {
 25     int radius = (diameter - 1) / 2;
 26     int imgW = srcImg.cols;
 27     int imgH = srcImg.rows;
 28     int channels = srcImg.channels();
 29     dstImg = Mat::zeros(imgH, imgW, CV_8UC1);
 30     int windowSize = diameter*diameter;
 31 
 32     //直方圖
 33     int Hist[256]={0};
 34     int histogramSize = 256;//灰度等級
 35     int thresholdValue = windowSize / 2 + 1;
 36 
 37     uchar *pSrcData=srcImg.data;
 38     uchar *pDstData=dstImg.data;
 39 
 40     int right=imgW-radius;
 41     int bot=imgH-radius;
 42     for (int i=radius; i<right; i++)
 43     {
 44         for (int j=radius; j<bot; j++)
 45         {
 46             //每一行第一個待濾波像素建立直方圖
 47             if(j==radius)
 48             {
 49                 memset(Hist, 0, histogramSize*sizeof(int));
 50                 for (int y=i-radius; y<=i+radius; y++)
 51                 {
 52                     for (int x=j-radius; x<=j+radius; x++)
 53                     {
 54                         uchar value=pSrcData[ y*imgW+x];
 55                         Hist[value]++;
 56                     }
 57                 }
 58             }
 59             else//更新直方圖
 60             {
 61                 int left=j-radius-1;
 62                 int right=j+radius;
 63                 for (int y=i-radius; y<=i+radius; y++)
 64                 {
 65                     //減去左邊一列
 66                     int leftIdx=y*imgW+left;
 67                     uchar leftValue=pSrcData[leftIdx];
 68                     Hist[leftValue]--;
 69 
 70                     //加上右邊一列
 71                     int rightIdx=y*imgW+right;
 72                     uchar rightValue=pSrcData[rightIdx];
 73                     Hist[rightValue]++;
 74                 }
 75             }
 76 
 77             //直方圖求中值
 78             uchar filterValue=calMediaValue(Hist, thresholdValue);
 79             pDstData[i*imgW+j]=filterValue;
 80         }
 81     }
 82 
 83     //邊界直接賦原始值，不做濾波處理
 84     pSrcData = srcImg.data;
 85     pDstData = dstImg.data;
 86     //上下邊界
 87     for (int j = 0; j < imgW; j++)
 88     {
 89         for (int i = 0; i < radius; i++)
 90         {
 91             int idxTop = i*imgW + j;
 92             pDstData[idxTop] = pSrcData[idxTop];
 93             int idxBot = (imgH - i - 1)*imgW + j;
 94             pDstData[idxBot] = pSrcData[idxBot];
 95         }
 96     }
 97     //左右邊界
 98     for (int i = radius; i < imgH - radius - 1; i++)
 99     {
100         for (int j = 0; j < radius; j++)
101         {
102             int idxLeft = i*imgW + j;
103             pDstData[idxLeft] = pSrcData[idxLeft];
104             int idxRight = i*imgW + imgW - j-1;
105             pDstData[idxRight] = pSrcData[idxRight];
106         }
107     }
108 }
109 
110 
111 void main()
112 {
113     string imgPath = "data/";
114     Mat srcImg = imread(imgPath + "moon.jpg", 0);
115     Mat dstImg;
116     double t0 = cv::getTickCount();
117     fastMedianBlur(srcImg, dstImg, 5);
118     double t1 = cv::getTickCount();
119     cout << "time=" << (t1 - t0) / cv::getTickFrequency() << endl;
120     
121     imwrite("data/test/srcImg.jpg", srcImg);
122     imwrite("data/test/myFilter.jpg", dstImg);
123 }

View Code

Code02:

  1 #include <opencv2\opencv.hpp>
  2 #include <iostream>
  3 #include <string>
  4 
  5 using namespace cv;
  6 using namespace std;
  7 
  8 //計算亮度中值和灰度<=中值的像素點個數
  9 void CalculateImage_MedianGray_PixelCount(const Mat &histogram, int pixelCount, int &medianValue, int &pixleCountLowerMedian)
 10 {
 11     float *data = (float *)histogram.data;//直方圖
 12     int sum = 0;
 13     for (int i = 0; i <= 255; ++i)
 14     {
 15         //
 16         sum += data[i];
 17         if (2 * sum>pixelCount)
 18         {
 19             medianValue = i;
 20             pixleCountLowerMedian = sum;
 21             break;
 22         }
 23     }
 24 }
 25 
 26 void fastMedianBlur(const Mat &srcImg, Mat &dstImg, int diameter)
 27 {
 28     int radius = (diameter - 1) / 2;
 29     int imgW = srcImg.cols;
 30     int imgH = srcImg.rows;
 31     int channels = srcImg.channels();
 32     dstImg = Mat::zeros(imgH, imgW, CV_8UC1);
 33     int windowSize = diameter*diameter;
 34     Mat windowImg = Mat::zeros(diameter, diameter, CV_8UC1);
 35 
 36     //直方圖
 37     Mat histogram;
 38     int histogramSize = 256;//灰度等級
 39     int thresholdValue = windowSize / 2 + 1;//step1.設置閾值(步驟參考：圖像的高效編程要點之四)
 40 
 41     //待處理圖像像素
 42     uchar *pDstData = dstImg.data + imgW*radius + radius;
 43     //整個圖像中窗口位置指針
 44     uchar *pSrcData = srcImg.data;
 45 
 46     //逐行遍歷
 47     for (int i = radius; i <= imgH - 1 - radius; i++)
 48     {
 49         //列指針
 50         uchar *pColDstData = pDstData;
 51         uchar *pColSrcData = pSrcData;
 52 
 53         //單個窗口指針
 54         uchar *pWindowData = windowImg.data;
 55         uchar *pRowSrcData = pColSrcData;
 56         //從源圖中提取窗口圖像
 57         for (int winy = 0; winy <= diameter - 1; winy++)
 58         {
 59             for (int winx = 0; winx <= diameter - 1; winx++)
 60             {
 61                 pWindowData[winx] = pRowSrcData[winx];
 62             }
 63             pRowSrcData += imgW;
 64             pWindowData += diameter;
 65         }
 66 
 67         //求直方圖,確定中值，並記錄亮度<=中值的像素點個數
 68         calcHist(&windowImg,
 69             1,//Mat的個數
 70             0,//用來計算直方圖的通道索引，通道索引依次排開
 71             Mat(),//Mat()返回一個空值，表示不用mask,
 72             histogram, //直方圖
 73             1, //直方圖的維數，如果計算2個直方圖，就為2
 74             &histogramSize, //直方圖的等級數(如灰度等級),也就是每列的行數
 75             0//分量的變化范圍
 76         );
 77 
 78         //求亮度中值和<=中值的像素點個數
 79         int medianValue, pixelCountLowerMedian;
 80         CalculateImage_MedianGray_PixelCount(histogram, windowSize, medianValue, pixelCountLowerMedian);
 81         ////////////滑動窗口操作結束///////////////////////
 82 
 83         //濾波
 84         pColDstData[0] = (uchar)medianValue;
 85 
 86         //處理同一行下一個像素
 87         pColDstData++;
 88         pColSrcData++;
 89         for (int j = radius + 1; j <= imgW - radius - 1; j++)
 90         {
 91             //維護滑動窗口直方圖
 92             //刪除左側
 93             uchar *pWinLeftData = pColSrcData - 1;
 94             float *pHistData = (float*)histogram.data;
 95             for (int winy = 0; winy < diameter; winy++)
 96             {
 97                 uchar grayValue = pWinLeftData[0];
 98                 pHistData[grayValue] -= 1.0;
 99                 if (grayValue <= medianValue)
100                 {
101                     pixelCountLowerMedian--;
102                 }
103                 pWinLeftData += imgW;
104             }
105 
106             //增加右側
107             uchar *pWinRightData = pColSrcData + diameter - 1;
108             for (int winy = 0; winy < diameter; winy++)
109             {
110                 uchar grayValue = pWinRightData[0];
111                 pHistData[grayValue] += 1.0;
112                 if (grayValue <= medianValue)
113                 {
114                     pixelCountLowerMedian++;
115                 }
116                 pWinRightData += imgW;
117             }
118             //計算新的中值
119             if (pixelCountLowerMedian > thresholdValue)
120             {
121                 while (1)
122                 {
123                     pixelCountLowerMedian -= pHistData[medianValue];
124                     medianValue--;
125                     if (pixelCountLowerMedian <= thresholdValue)
126                     {
127                         break;    
128                     }
129                 }
130             }
131             else 
132             {
133                 while (pixelCountLowerMedian < thresholdValue)
134                 {
135                     medianValue++;
136                     pixelCountLowerMedian += pHistData[medianValue];
137 
138                 }
139             }
140 
141             pColDstData[0] = medianValue;
142             //下一個像素
143             pColDstData++;
144             pColSrcData++;
145         }
146         //移動至下一行
147         pDstData += imgW;
148         pSrcData += imgW;
149     }
150 
151     //邊界直接賦原始值，不做濾波處理
152     pSrcData = srcImg.data;
153     pDstData = dstImg.data;
154     //上下邊界
155     for (int j = 0; j < imgW; j++)
156     {
157         for (int i = 0; i < radius; i++)
158         {
159             int idxTop = i*imgW + j;
160             pDstData[idxTop] = pSrcData[idxTop];
161             int idxBot = (imgH - i - 1)*imgW + j;
162             pDstData[idxBot] = pSrcData[idxBot];
163         }
164     }
165     //左右邊界
166     for (int i = radius; i < imgH - radius - 1; i++)
167     {
168         for (int j = 0; j < radius; j++)
169         {
170             int idxLeft = i*imgW + j;
171             pDstData[idxLeft] = pSrcData[idxLeft];
172             int idxRight = i*imgW + imgW - j-1;
173             pDstData[idxRight] = pSrcData[idxRight];
174         }
175     }
176 }
177 
178 
179 void main()
180 {
181     string imgPath = "data/";
182     Mat srcImg = imread(imgPath + "moon.jpg", 0);
183     Mat dstImg;
184     double t0 = cv::getTickCount();
185     fastMedianBlur(srcImg, dstImg, 5);
186     //cv::medianBlur(srcImg, dstImg, 5); //OpenCV
187     double t1 = cv::getTickCount();
188     cout << "time=" << (t1 - t0) / cv::getTickFrequency() << endl;
189     
190     imwrite("data/test/srcImg.bmp", srcImg);
191     imwrite("data/test/myFilter.bmp", dstImg);
192 }

View Code

免責聲明！

本站轉載的文章為個人學習借鑒使用，本站對版權不負任何法律責任。如果侵犯了您的隱私權益，請聯系本站郵箱yoyou2525@163.com刪除。

猜您在找 verilog 實現中值濾波中值濾波 java實現 verilog實現中值濾波基於FPGA的中值濾波算法實現基於MATLAB的中值濾波算法實現中值濾波的幾種快速算法基於MATLAB的中值濾波均值濾波以及高斯濾波的實現中值濾波中值濾波 FPGA實現圖像的非線性濾波：中值濾波