《圖像處理實例》之尋找山脊線

時間 2019-11-16

標籤圖像處理實例尋找脊線简体版

原文原文鏈接

尋找山脊線算法ios

如下的改進是http://www.imagepy.org/的做者原創，我只是對其理解以後改進和說明，歡迎你們使用這個小軟件！
c++

C++版本已經編寫完畢，有須要的朋友請留郵箱算法

算法含義：函數

　　　　提取一個山脈的山峯，圖像中就是在距離變換以後提取局部最亮的線ui

算法應用：spa

　　　　用做圖像的斷裂修補3d

方法：code

　　　　1.直接用分水嶺算法進行區域分割。orm

　　　　2.相似分水嶺的算法，找山脊線。對象

分水嶺算法：

　　　　這裏有個哲學的思想，仍是開頭大佬點撥個人，修補前景就等價於分割背景，必定要理解這句話的含義，咱們把白色的直線當作前景修復對象，那麼就等價於用白線去分割黑色的背景。

　　　　opencv的分水嶺算法不是很好，必須找到合適的mask才能達到要求。

　　　　skimage的分水嶺相對較好，效果差強人意吧.

　　　　　　分水嶺結果如上圖，在這基礎上閾值一下就能夠獲得咱們的圖像（這個圖像有多餘的線條），而後再和原始圖像作差值獲得修補的線條，把修補的很長線條去除以後再加上原圖就能夠了。

　　　　你們也發現了一個問題，手動的元素比較多，預處理作的很差（局部極大值找的很差，我用的自適應閾值，效果很差就不上圖了），上圖的效果算比較好的，至少能夠解決咱們的問題。

代碼：

  1 #include<opencv2/opencv.hpp>
  2 #include<iostream>
  3 using namespace std;
  4 using namespace cv;
  5 
  6 #define WINDOW_NAME "【程序窗口1】"
  7 Mat g_maskImage, g_srcImage;
  8 Point prevPt(-1, -1);
  9 static void on_Mouse(int event, int x, int y, int flags, void*);
 10 
 11 int main()
 12 {
 13     //載入原圖，初始化掩膜和灰度圖
 14     g_srcImage = imread("123.png");
 15     imshow(WINDOW_NAME, g_srcImage);
 16     Mat srcImage, grayImage;
 17     g_srcImage.copyTo(srcImage);
 18     cvtColor(g_srcImage, g_maskImage, COLOR_BGR2GRAY);
 19     cvtColor(g_maskImage, grayImage, COLOR_GRAY2BGR);
 20     g_maskImage = Scalar::all(0);
 21     //設置鼠標回調函數
 22     setMouseCallback(WINDOW_NAME, on_Mouse);
 23     //輪詢按鍵
 24     while (1)
 25     {
 26         int c = waitKey(0);
 27         if ((char)c == 27) break;
 28         if ((char)c == '2') {   //按鍵‘2’， 恢復源圖
 29             g_maskImage = Scalar::all(0);
 30             srcImage.copyTo(g_srcImage);
 31             imshow("image", g_srcImage);
 32         }
 33         if ((char)c == '1' || (char)c == ' ') {
 34             //定義一些參數            
 35             vector<vector<Point> > contours;
 36             vector<Vec4i> hierarchy;
 37             //尋找輪廓
 38             findContours(g_maskImage, contours, hierarchy, CV_RETR_CCOMP, CV_CHAIN_APPROX_SIMPLE);
 39             //輪廓爲空時的處理
 40             if (contours.empty()) continue;
 41             //複製掩膜
 42             Mat maskImage(g_maskImage.size(), CV_32S);
 43             maskImage = Scalar::all(0);
 44             //循環繪製出輪廓
 45             int compCount = 0;
 46             for (int index = 0; index >= 0; index = hierarchy[index][0], compCount++)
 47                 drawContours(maskImage, contours, index, Scalar::all(compCount + 1), -1, LINE_8, hierarchy);
 48             //compCount爲零時的處理
 49             if (compCount == 0)
 50                 continue;
 51             //生成隨機顏色
 52             vector<Vec3b> colorTab;
 53             for (unsigned int i = 0; i < compCount; i++) {
 54                 int b = theRNG().uniform(0, 255);
 55                 int g = theRNG().uniform(0, 255);
 56                 int r = theRNG().uniform(0, 255);
 57                 colorTab.push_back(Vec3b((uchar)b, (uchar)g, (uchar)r));
 58             }
 59             //計算處理時間並輸出到窗口中
 60             double dTime = (double)getTickCount();
 61             watershed(srcImage, maskImage);
 62             dTime = (double)getTickCount() - dTime;
 63             printf("處理時間=%gms\n", dTime*1000. / getTickFrequency());
 64             //雙層循環，將分水嶺圖像遍歷存入watershedImage中
 65             Mat watershedImage(maskImage.size(), CV_8UC3);
 66             for (unsigned int i = 0; i < maskImage.rows; i++)
 67                 for (unsigned int j = 0; j < maskImage.cols; j++)
 68                 {
 69                     int index = maskImage.at<int>(i, j);
 70                     if (index == -1)
 71                         watershedImage.at<Vec3b>(i, j) = Vec3b(255, 255, 255);
 72                     else if (index <= 0 || index > compCount)
 73                         watershedImage.at<Vec3b>(i, j) = Vec3b(0, 0, 0);
 74                     else
 75                         watershedImage.at<Vec3b>(i, j) = colorTab[index - 1];
 76                 }
 77             //混合灰度圖和分水嶺效果圖並顯示最終的窗口
 78             //watershedImage = watershedImage*0.5 + grayImage*0.5;
 79             imshow("watershed transform", watershedImage);
 80         }
 81     }
 82     return 0;
 83 }
 84 
 85 static void on_Mouse(int event, int x, int y, int flags, void*)
 86 {
 87     //處理鼠標不在窗口中的狀況
 88     if (x < 0 || x >= g_srcImage.cols || y < 0 || y >= g_srcImage.rows) return;
 89     //處理鼠標左鍵相關消息
 90     if (event == EVENT_LBUTTONUP || !(flags & EVENT_FLAG_LBUTTON))   //左鍵擡起動做或處於沒有按下狀態
 91         prevPt = Point(-1, -1);
 92     else if (event == EVENT_LBUTTONDOWN)  //左鍵按下動做
 93         prevPt = Point(x, y);
 94     //鼠標左鍵按下並移動，繪製出白色線條
 95     else if (event == EVENT_MOUSEMOVE && (flags & EVENT_FLAG_LBUTTON))
 96     {
 97         Point pt(x, y);
 98         if (prevPt.x < 0) prevPt = pt;
 99         line(g_maskImage, prevPt, pt, Scalar::all(255), 5);
100         line(g_srcImage, prevPt, pt, Scalar::all(255), 5);
101         prevPt = pt;
102         imshow(WINDOW_NAME, g_srcImage);
103     }
104 
105 }

找山脊線算法：

　　　　　相似於骨架提取、圖像細化、中軸線提取等思想

　　　把白色的線當作前景圖，而後進行距離變換，以後進行全局閾值把前景所有包括在mask裏面，接着利用原始圖對mask進行操做，具體的細節看代碼分析。

　　　　如下是效果圖，很完美的解決了問題。

算法的原理：

　　　　第一步：對圖像進行距離變換---->>>以後閾值化（全局閾值就好）

例如：效果圖

　　　　註釋：這裏不進行距離變換，直接進行閾值化也能夠，可是對於特殊的圖像效果不好。由於掩膜的應用能夠結合分水嶺[0:255]進行操做，沒有距離變換的掩膜就只有255，這就和圖像細化或者中軸線算法差很少了。

原圖

距離變換圖

掩膜圖

結果圖

掩膜圖

結果圖

　　　　第二步：對二值化圖像進行標記（作mask）

　　　　　　　　邊界標記成4，這裏爲了防止越界，按照opencv c++的說法就是img.rows-1....

　　　　　　　　背景標記成2

　　　　　　　　前景標記成0，這裏是除去圖形邊界的前景

　　　　　　　　圖形邊界標記成1

　　　　第三步：分水嶺算法

　　　　　　　　進行操做的是距離圖+掩膜圖（dis+mask）

　　　　　　　　這裏使用漲水的方式對像素進行操做，從剛開始標記的邊界1開始向內部漲水，因爲距離變換以後映射到【0-255】，因此水從0開始漲，知道255爲止：

　　　　　　　　這一步關鍵一點是內部的凹池，採用相似水漫算法進行：

　　　　　　　　　　就這樣把低於水平線的像素都填充爲背景2了，固然按照其餘的行列循環也是能夠的，怎麼寫都行！

　　　　註釋：

　　　　　　這裏有幾個小技巧：

　　　　　　　　　　　　1.由於是從邊界開始判斷的，好比邊界高度爲100，那麼水就從100開始漲，否則0-99是沒意義的，具體代碼見程序說217行。

　　　　　　　　　　　　2.每一次的循環最少是圖形一圈，由於這裏以漲水爲標準，好比上面的以原圖二值化爲掩膜，那麼圖像都爲255，這時候就是一次性得出山脊線。

　　　　　　　　　　　　　簡單的說就是每次處理的單位都是一個水平線，好比水平線是200，那麼邊界高度爲201就是下一次漲水再進行處理。

　　　　　　　　　　　　3.掩膜取得合適得出的圖像效果最好，也不是越大越好，具體例子就不舉了。

　　　　　　　　　　　　4.由於這裏查表遵循的是不斷裂原則（貌似圖像細化和骨架提取都是這樣），回頭想一想你把掩膜圖進行細化而後再和原圖疊加也能進行圖像修復哈！

　　　　下面是我本身的一個誤區：（這裏你們看完上面說明就知道如下的問題是不存在的）　　　　

　　　　　　　　固然這裏的算法不穩定，一些特殊的圖像仍是不行，好比下面的圖像，水還沒漲到200，圖就已經分割完了，程序是水漲一次就分割一次，可是這個圖明顯都大於200.

  1 from scipy.misc import imread
  2 import matplotlib.pyplot as plt
  3 import numpy as np
  4 from numba import jit
  5 from skimage.data import camera
  6 import scipy.ndimage as ndimg
  7 from time import time
  8 import cv2
  9 from scipy.ndimage import label, generate_binary_structure
 10 
 11 strc = np.ones((3, 3), dtype=np.bool)
 12 
 13 
 14 def count(n):
 15     a = [(n>>i) & 1 for i in range(8)]
 16     if sum(a)<=1:return False
 17     if a[1] & a[3] & a[5] & a[7]:return False
 18     a = np.array([[a[0],a[1],a[2]],
 19                   [a[7],  0 ,a[3]],
 20                   [a[6],a[5],a[4]]])
 21     n = label(a, strc)[1]
 22     return n<2
 23 
 24 
 25 lut = np.array([count(n) for n in range(256)])
 26 lut = np.dot(lut.reshape((-1, 8)), [1, 2, 4, 8, 16, 32, 64, 128]).astype(np.uint8)
 27 
 28 
 29 #@jit
 30 def core(n):
 31     a = np.zeros(8, dtype=np.uint8)
 32     for i in range(8):
 33         a[i] = (n >> i * 2) & 3
 34     if a[1] == 1 and a[0] == 0: a[0] = 1
 35     if a[1] == 1 and a[2] == 0: a[0] = 1
 36     if a[3] == 1 and a[2] == 0: a[2] = 1
 37     if a[3] == 1 and a[4] == 0: a[4] = 1
 38     if a[5] == 1 and a[4] == 0: a[4] = 1
 39     if a[5] == 1 and a[6] == 0: a[6] = 1
 40     if a[7] == 1 and a[6] == 0: a[6] = 1
 41     if a[7] == 1 and a[0] == 0: a[0] = 1
 42     for i in range(8):
 43         if a[i] == 0 or a[i] == 2: a[i] = 0
 44         if a[i] == 1 or a[i] == 3: a[i] = 1
 45     return np.dot(a, [1, 2, 4, 8, 16, 32, 64, 128])
 46 
 47 
 48 index = np.array([core(i) for i in range(65536)], dtype=np.uint8)
 49 
 50 #
 51 '''
 52 lut = np.array([223, 221, 1, 221, 1, 221, 1, 221, 1, 0, 0, 0, 1, 221, 1, 221, 207, 204,
 53                 0, 204, 207,  51, 207, 1, 207, 204, 0, 204, 207, 51, 207, 51], dtype=np.uint8)
 54 '''
 55 
 56 #
 57 def nbs8(h, w):
 58     return np.array([-w - 1, -w, -w + 1, +1, +w + 1, +w, +w - 1, -1], dtype=np.int32)
 59 
 60 #相似建立2X2內核（上下左右4個數）
 61 def nbs4(h, w):
 62     return np.array([-1, -w, 1, w], dtype=np.int32)
 63 
 64 
 65 #@jit
 66 
 67 def fill(img, msk, p, level, pts, s, nbs, buf):
 68     n = 0
 69     cur = 0
 70     buf[0] = p      
 71     msk[p] = 2      
 72     bs = 1
 73     while cur < bs:
 74         p = buf[cur]
 75         for dp in nbs:
 76             if msk[p + dp] != 0: continue  
 77             if img[p + dp] < level: 
 78                 buf[bs] = p + dp            
 79                 msk[p + dp] = 2             
 80                 bs += 1                     
 81                 if bs == len(buf):
 82                     buf[:bs - cur] = buf[cur:bs]
 83                     bs -= cur
 84                     cur = 0
 85             else:                           
 86                 pts[s + n] = p + dp
 87                 msk[ p + dp] = 1
 88                 n += 1
 89         cur += 1
 90     return n
 91 
 92 
 93 #@jit
 94 def check(msk, p, nbs, lut):
 95     c = 0
 96     s = 0
 97     for i in range(8):
 98         v = msk[p + nbs[i]]
 99         # if v==0: c|=(0<<i*2)
100         if v == 1: c |= (1 << i * 2)
101         if v == 2: c |= (2 << i * 2)
102         if v == 3: c |= (3 << i * 2)
103     v = index[c]
104     if lut[v // 8] >> v % 8 & 1:
105         msk[p] = 2
106     else:
107         msk[p] = 3
108 
109 
110 #@jit
111 def step(img, msk, pts, s, level, nbs, nbs8):
112     ddd = 0
113     cur = 0
114     buf = np.zeros(10240, dtype=np.int64)
115     while cur < s:  
116         p = pts[cur]    
117         if img[p] > level:
118             cur += 1
119             continue
120 
121         filled = False
122         for dp in nbs:
123             if msk[p + dp] == 4: msk[p] = 2 
124             if msk[p + dp] == 0:            
125                 if img[p + dp] >= level:    
126                     msk[p + dp] = 1         
127                     pts[s] = p + dp         
128                     s += 1                  
129                 else:                       
130                     n = fill(img, msk, p, level, pts, s, nbs, buf)
131                     s += n      
132                     filled = True
133 
134         if filled:      
135             cur += 1
136             continue
137         elif msk[p] == 1:  
138             check(msk, p, nbs8, lut)   
139 
140         cur += 1   
141     return cur
142 
143 
144 @jit
145 def clear(msk, pts, s):
146     cur = 0
147     for c in range(s):
148         if msk[pts[c]] == 1:
149             pts[cur] = pts[c]
150             cur += 1
151     return cur
152 
153 
154 @jit
155 def collect(img, mark, nbs, pts):
156     bins = np.zeros(img.max() + 1, dtype=np.uint32) 
157     cur = 0     
158     '''
159     for p in range(len(mark)):      
160         if mark[p] != 0: continue
161         for dp in nbs:          
162             if mark[p + dp] == 1:
163                 mark[p] = 2
164     '''
165     for p in range(len(mark)):         
166         if mark[p] == 1: mark[p] = 2
167     for p in range(len(mark)):
168         if mark[p] != 0: continue
169         s = 0   
170         for dp in nbs:
171             if mark[p + dp] == 2:       
172                 s += 1                  
173                 mark[p] = 1             
174                 pts[cur] = p           
175                 cur += 1               
176                 break
177         if s == 0: bins[img[p]] += 1    
178     return cur, bins
179 
180 
181 #@jit
182 def watershed(img, mark):
183     oimg, omark = img, mark
184     ndim = img.ndim
185     mark[[0, -1], :] = 4   
186     mark[:, [0, -1]] = 4    
187 
188     nb4 = nbs4(*img.shape) 
189     nb8 = nbs8(*img.shape)
190     acc = np.cumprod((1,) + img.shape[::-1][:-1])[::-1] 
191     img = img.ravel()       
192     mark = mark.ravel()     
193 
194     pts = np.zeros(131072, dtype=np.int64) 
195     s, bins = collect(img, mark, nb4, pts)  
196     for level in range(len(bins)):  
197         if bins[level] == 0: continue
198         s = clear(mark, pts, s)     
199         s = step(img, mark, pts, s, level, nb4, nb8)
200 '''
201 dem = cv2.imread('gis1.jpg',0)
202 plt.imshow(dem, cmap='gray')
203 mark = (dem > 150).astype(np.uint8)
204 plt.imshow(mark, cmap='gray')
205 plt.show()
206 
207 watershed(dem, mark.copy())
208 start = time()
209 watershed(dem, mark)
210 print(time() - start)
211 dem[:] = dem * 0.8
212 dem[mark == 3] = 255
213 plt.imshow(dem, cmap='gray')
214 plt.show()
215 '''
216 #dem = imread('line2.png')
217 dem = cv2.imread('0.jpg')
218 dem = cv2.cvtColor(dem,cv2.COLOR_BGR2GRAY)
219 ret2, dem = cv2.threshold(dem, 0, 255, cv2.THRESH_BINARY | cv2.THRESH_OTSU)
220 dis = ndimg.distance_transform_edt(~dem).astype(np.uint8)
221 dis = ~dis
222 mark = (dis<230).astype(np.uint8)
223 plt.imshow(dis, cmap='gray')
224 plt.imshow(mark, cmap='gray')
225 watershed(dis, mark)
226 dem //= 2
227 dem[mark==3] = 255
228 plt.imshow(dem, cmap='gray')
229 plt.show()

C++版本：

  1 //---table
  2 uchar lut[] = { 200, 206, 220, 204, 0, 207, 0, 204, 0, 207, 221, 51, 1, 207, 221, 51,
  3 0, 0, 221, 204, 0, 0, 0, 204, 1, 207, 221, 51, 1, 207, 221, 51 };
  4 
  5 //---FindRidge
  6 void FindRidge(InputArray& _src, Mat& mask, vector<Point> Edge_Point, uchar thred)
  7 {
  8     Mat src = _src.getMat();// , mask = _mask.getMat();
  9 
 10     mask = src.clone();
 11     bitwise_not(mask, mask);
 12     distanceTransform(mask, mask, DIST_L2, DIST_MASK_3, 5);
 13     normalize(mask, mask, 0, 255, NORM_MINMAX);
 14     mask.convertTo(mask, CV_8UC1);
 15     threshold(mask, mask, thred, 255, THRESH_BINARY_INV);
 16     Mat temp = mask.clone();
 17 
 18     bitwise_not(mask, mask);
 19     mask = mask - 253;
 20     Mat rows = Mat::ones(Size(src.cols, 1), CV_8UC1), cols = Mat::zeros(Size(1, src.rows), CV_8UC1);
 21     rows.setTo(4); cols.setTo(4);
 22     Mat src_rows_beg = mask.row(0),                src_cols_beg = mask.col(0);
 23     Mat src_rows_end = mask.row(src.rows - 1),    src_cols_end = mask.col(src.cols - 1);
 24     rows.copyTo(src_rows_beg); rows.copyTo(src_rows_end);
 25     cols.copyTo(src_cols_beg); cols.copyTo(src_cols_end);
 26 
 27     Edge_Point = FindEdge(temp, mask);
 28     bool TF = true;
 29     while (TF)
 30     {
 31         TF = false;
 32         for (size_t i = 0; i < mask.rows - 1; i++)
 33         {
 34             uchar* msk_up = mask.ptr(i - 1);
 35             uchar* msk = mask.ptr(i);
 36             uchar* msk_dw = mask.ptr(i + 1);
 37             for (size_t j = 0; j < mask.cols - 1; j++)
 38             {
 39                 uchar _temp_data = msk[j];
 40                 msk[j] = msk[j] == 2 && (msk_up[j] == 255 || msk_up[j] == 0)
 41                     && (msk_dw[j] == 255 || msk_dw[j] == 0)
 42                     && (msk[j - 1] == 255 || msk[j - 1] == 0)
 43                     && (msk[j + 1] == 255 || msk[j + 1] == 0) ? 0 : msk[j];
 44                 msk[j] = msk[j] == 0 && (msk_up[j] == 2 || msk_dw[j] == 2 || msk[j - 1] == 2 || msk[j + 1] == 2) ? 2 : msk[j];
 45                 TF = _temp_data != msk[j] ? true : TF;
 46             }
 47         }
 48     }
 49 
 50     distanceTransform(src, src, DIST_L2, DIST_MASK_3, 5);
 51     normalize(src, src, 0, 255, NORM_MINMAX);
 52     src.convertTo(src, CV_8UC1);
 53 
 54     const int histSize = 255;
 55     float range[] = { 0, 255 };
 56     const float* histRange = { range };
 57     Mat hist;
 58     calcHist(&src, 1, 0, Mat(), hist, 1, &histSize, &histRange, true, false);
 59 
 60     hist = hist.reshape(1,1);
 61     normalize(hist, hist, 0, 1000, NORM_MINMAX);
 62     hist.convertTo(hist, CV_32FC1);
 63     //float* ptr = hist.ptr(0);
 64     for (size_t level = 0; level <= 255; level++)
 65     {
 66         if (!hist.at<float>(0,level)) continue;
 67         FloorEdge(src, Edge_Point, mask, level);
 68     }
 69 }
 70 
 71 /*
 72 
 73 vector<Point> FindEdge(InputArray& _src, Mat& mask)
 74 {
 75     Mat src1 = _src.getMat(),src = src1.clone();
 76     Mat kernel = getStructuringElement(MORPH_RECT, Size(3, 3));
 77     morphologyEx(src, src, MORPH_ERODE, kernel);
 78     vector<Point> wjy_Point;
 79     vector<vector<Point>> contours;
 80     vector<Vec4i> hierarchy;
 81     findContours(src, contours, hierarchy, RETR_TREE, CHAIN_APPROX_NONE, Point(-1, -1));
 82     for (size_t i = 0; i < contours.size(); i++)
 83     {
 84         for (size_t j = 0; j < contours[i].size(); j++)
 85         {
 86             wjy_Point.push_back(contours[i][j]);
 87             mask.at<uchar>(contours[i][j]) = 255;
 88         }
 89     }
 90     return wjy_Point;
 91 }
 92 */
 93 
 94 
 95 vector<Point> FindEdge(Mat& src,Mat& mask)
 96 {
 97     //Mat src = mask.clone();
 98     Mat kernel = getStructuringElement(MORPH_RECT, Size(3, 3));
 99     //morphologyEx(src, src, MORPH_ERODE, kernel);
100     vector<Point> wjy_Point;
101     vector<vector<Point>> contours;
102     vector<Vec4i> hierarchy;
103     findContours(src, contours, hierarchy, RETR_TREE, CHAIN_APPROX_NONE, Point(-1, -1));
104     for (size_t i = 0; i < contours.size(); i++)
105     {
106         for (size_t j = 0; j < contours[i].size(); j++)
107         {
108             wjy_Point.push_back(contours[i][j]);
109             mask.at<uchar>(contours[i][j]) = 255;
110         }
111     }
112     return wjy_Point;
113 }
114 
115 
116 void FloorEdge(InputArray& _src, vector<Point>& Edge_Point, Mat& mask,int level)
117 {
118     Mat src = _src.getMat();
119     for (int i = 0; i < Edge_Point.size(); i++)
120     {
121         vector<Point> temp_vector;
122         temp_vector.push_back(Point(Edge_Point[i].x, Edge_Point[i].y - 1));
123         temp_vector.push_back(Point(Edge_Point[i].x, Edge_Point[i].y + 1));
124         temp_vector.push_back(Point(Edge_Point[i].x - 1, Edge_Point[i].y));
125         temp_vector.push_back(Point(Edge_Point[i].x + 1, Edge_Point[i].y));
126         uchar* msk = mask.ptr(Edge_Point[i].y);
127         uchar* img = src.ptr(Edge_Point[i].y);
128         if (img[Edge_Point[i].x] > level)    continue;
129         bool Flag = true;
130         uchar count_num = 0;
131         for (size_t j = 0; j < temp_vector.size(); j++)
132         {
133             uchar* pre_data = mask.ptr(temp_vector[j].y);
134             if (pre_data[temp_vector[j].x] == 2 || pre_data[temp_vector[j].x] == 4)
135             {
136                 pre_data[temp_vector[j].x] = 2;
137                 continue;
138             }
139             else if (pre_data[temp_vector[j].x] == 128 || pre_data[temp_vector[j].x] == 255)
140             {
141                 continue;                
142                 count_num++;
143             }
144             else
145             {
146                 if (src.at<uchar>(temp_vector[j]) >= level)
147                 {            
148                     count_num++;
149                     pre_data[temp_vector[j].x] = 255;
150                     Edge_Point.push_back(temp_vector[j]);    
151                     //Edge_Point.insert(Edge_Point.begin()+i+1,temp_vector[j]);
152                 }
153                 else
154                 {
155                     int temp = FillBlock(src, Edge_Point, mask, level, Edge_Point[i]);
156                     Flag = false;
157                 }
158             }
159         }
160         if (false)
161         {
162             msk[Edge_Point[i].x] = 128;
163             Edge_Point.erase(Edge_Point.begin() + i);
164             i--;
165             continue;
166         }
167         else if (count_num == 4)
168         {
169             Edge_Point.push_back(Edge_Point[i]);
170             Edge_Point.erase(Edge_Point.begin() + i);
171             i--;
172             continue;
173         }
174         else if (msk[Edge_Point[i].x] == 255)
175         {
176             msk[Edge_Point[i].x] = CheckData(mask, Edge_Point[i]) == true ? 128 : 2;
177             Edge_Point.erase(Edge_Point.begin() + i);
178             i--;
179         }
180     }
181 }
182 
183 int FillBlock(InputArray& _src, vector<Point>& Edge_Point, Mat& mask, int level, Point center)
184 {
185     Mat src = _src.getMat();
186     mask.at<uchar>(center) = 2;
187     vector<Point> fill_point;
188     int count = 0, count_mount = 1;
189     fill_point.push_back(center);
190     while (count < count_mount)
191     {
192         vector<uchar*> img;
193         vector<uchar*> msk;
194         for (int i = -1; i < 2; i++)
195         {
196             img.push_back(src.ptr<uchar>(fill_point[count].y + i));
197             msk.push_back(mask.ptr<uchar>(fill_point[count].y + i));
198         }
199         for (size_t i = 0; i < 3; i++)
200         {
201             for (int j = -1; j < 2; j++)
202             {
203                 if (img[i][fill_point[count].x + j] < level && !(j == 0 && i == 1) && msk[i][fill_point[count].x + j] == 0)
204                 {
205                     fill_point.push_back(Point(fill_point[count].x + j, fill_point[count].y + i - 1));
206                     msk[i][fill_point[count].x + j] = 2;
207                 }
208                 else if (img[i][fill_point[count].x + j] >= level && msk[i][fill_point[count].x + j] == 0)
209                 {
210                     Edge_Point.push_back(Point(fill_point[count].x + j, fill_point[count].y + i - 1));
211                     msk[i][fill_point[count].x + j] = 255;
212                 }
213             }
214         }
215         //msk[1][fill_point[count].x] = 2;
216         count_mount = fill_point.size() - 1;
217         fill_point.erase(fill_point.begin());
218     }
219     return 0;
220 }
221 
222 
223 bool CheckData(Mat& mask, Point center)
224 {
225     uchar* msk_up = mask.ptr(center.y - 1);
226     uchar* msk = mask.ptr(center.y);
227     uchar* msk_dw = mask.ptr(center.y + 1);
228     int num[8];
229     int sum = (num[0] = msk_up[center.x - 1] == 255 || msk_up[center.x - 1] == 128 || msk_up[center.x - 1] == 0 ? 1 : 0)
230             + (num[1] = msk_up[center.x] == 255 || msk_up[center.x] == 128 || msk_up[center.x] == 0 ? 1 : 0) * 2
231             + (num[2] = msk_up[center.x + 1] == 255 || msk_up[center.x + 1] == 128 || msk_up[center.x + 1] == 0 ? 1 : 0) * 4
232             + (num[3] = msk[center.x - 1] == 255 || msk[center.x - 1] == 128 || msk[center.x - 1] == 0 ? 1 : 0) * 8
233             + (num[4] = msk[center.x + 1] == 255 || msk[center.x + 1] == 128 || msk[center.x + 1] == 0 ? 1 : 0) * 16
234             + (num[5] = msk_dw[center.x - 1] == 255 || msk_dw[center.x - 1] == 128 || msk_dw[center.x - 1] == 0 ? 1 : 0) * 32
235             + (num[6] = msk_dw[center.x] == 255 || msk_dw[center.x] == 128 || msk_dw[center.x] == 0 ? 1 : 0) * 64
236             + (num[7] = msk_dw[center.x + 1] == 255 || msk_dw[center.x + 1] == 128 || msk_dw[center.x + 1] == 0 ? 1 : 0) * 128;
237     return ((lut[uchar(sum / 8)] >> sum % 8) & 1) != 1 ? true : false;
238 }

相關標籤/搜索

每日一句

每一个你不满意的现在，都有一个你没有努力的曾经。

《圖像處理實例》 之 尋找山脊線

《圖像處理實例》之尋找山脊線