1 // Copyright (C) 2010 Andrew Kelly, IPS Radio & Space Services 2 // 3 // This program is free software; you can redistribute it and/or modify it under 4 // the terms of the GNU General Public License as published by the Free Software 5 // Foundation; either version 3 of the License, or (at your option) any later 6 // version. 7 // 8 // This program is distributed in the hope that it will be useful, but WITHOUT 9 // ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or 10 // FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more 11 // details. 12 // 13 // You should have received a copy of the GNU General Public License along with 14 // this program; if not, see <http://www.gnu.org/licenses/>. 15 16 #include <vector> 17 18 /** 19 * Oct-file to trace the boundary of an object in a binary image. 20 * 21 * b = boundary(region, conn=8) 22 */ 23 #include <octave/oct.h> 24 25 using namespace std; 26 27 DEFUN_DLD(__boundary__, args, nargout, 28 "-*- texinfo -*-\n\ 29 @deftypefn {Loadable Function} {} boundary(@var{region})\n\ 30 @deftypefnx {Loadable Function} {} boundary(@var{region}, @var{conn})\n\ 31 Trace the boundary of an object in a binary image.\n\ 32 \n\ 33 @code{boundary} computes the exterior clockwise boundary of the single \ 34 @var{conn}-connected object represented by the non-zero pixels \ 35 of @var{region}. It uses an algorithm based on Moore-neighbour tracing.\n\ 36 \n\ 37 @var{conn} can be either 8 (the default) or 4.\n\ 38 \n\ 39 @var{b} is an N-by-2 matrix containing the row/column coordinates of points \ 40 on the boundary. The first boundary point is the first non-zero \ 41 pixel of @var{region}, as determined by @code{find}. The last boundary \ 42 point is the same as the first.\n\ 43 @seealso{boundaries, bwlabel, find}\n\ 44 @end deftypefn") 45 { 46 octave_value_list retval; 47 48 enum { ROW, COL }; 49 50 // check number of arguments 51 const int nargin = args.length (); 52 if (nargin > 2 || nargout != 1) 53 error ("__boundary__: wrong number of input arguments"); 54 55 // extract arguments 56 const boolMatrix unpadded = args (0).bool_matrix_value (); 57 const int conn = (nargin > 1) 58 ? (int) args (1).scalar_value () 59 : 8; 60 61 // pad to avoid boundary issues 62 int rows = unpadded.rows (); 63 int cols = unpadded.columns (); 64 boolMatrix region (rows + 2, cols + 2, false); 65 for (int r = 0; r < rows; r++) 66 for (int c = 0; c < cols; c++) 67 region.elem (r+1, c+1) = unpadded (r, c); 68 69 // the padded size 70 rows += 2; 71 cols += 2; 72 73 // find the (first two) true pixels, if any 74 std::vector <int> pixels; 75 for (int i = 0; pixels.size () < 2 && i < region.numel (); ++i) 76 if (region.elem (i)) 77 pixels.push_back (i); 78 79 if (pixels.empty ()) 80 return retval; 81 82 // the starting boundary point 83 const int start = pixels [0]; 84 std::vector <int> bound; 85 bound.push_back (start); 86 87 // is this the only point? 88 if (pixels.size () == 1) 89 bound.push_back (start); 90 91 // otherwise, find the boundary by tracing the Moore neighbourhood of its pixels 92 // 93 // 8-connected: 7 0 1 4-connected: 0 94 // 6 . 2 3 . 1 95 // 5 4 3 2 96 else 97 { 98 // relative row/column positions 99 static const int row8 [] = {-1, -1, 0, 1, 1, 1, 0, -1}; 100 static const int col8 [] = { 0, 1, 1, 1, 0, -1, -1, -1}; 101 static const int row4 [] = {-1, 0, 1, 0 }; 102 static const int col4 [] = { 0, 1, 0, -1 }; 103 const int* mr = (conn == 4) ? row4 : row8; 104 const int* mc = (conn == 4) ? col4 : col8; 105 106 // next after backing-up 107 static const int back8 [] = {7, 7, 1, 1, 3, 3, 5, 5}; 108 static const int back4 [] = {3, 0, 1, 2}; 109 const int* mBack = (conn == 4) ? back4 : back8; 110 111 // relative indexes into the region for the Moore neighbourhood pixels 112 std::vector<int> mi (conn); 113 for (int i = 0; i < conn; ++i) 114 mi[i] = mr[i] + (rows * mc [i]); 115 116 // next neighbourhood pixel 117 static const int next8 [] = {1, 2, 3, 4, 5, 6, 7, 0}; 118 static const int next4 [] = {1, 2, 3, 0}; 119 const int* mNext = (conn == 4) ? next4 : next8; 120 121 // the final boundary point to be visited 122 int finish = 0; 123 for (int i = 0; i < conn; ++i) 124 if (region.elem(start + mi[i])) 125 finish = start + mi[i]; 126 127 // look for the next boundary point, starting at the next neighbour 128 int bp = start; 129 int mCurrent = mNext [0]; 130 bool done = false; 131 while (!done) 132 { 133 // next neighbour 134 int cp = bp + mi[mCurrent]; 135 136 // if this pixel is false, try the next one 137 if (!region.elem (cp)) 138 { 139 mCurrent = mNext [mCurrent]; 140 } 141 // otherwise, we have another boundary point 142 else 143 { 144 bound.push_back (cp); 145 146 // either we're back at the start for the last time 147 if (bp == finish && cp == start) 148 { 149 done = true; 150 } 151 // or we step back to where we came in from, and continue 152 else 153 { 154 bp = cp; 155 mCurrent = mBack [mCurrent]; 156 } 157 } 158 } 159 } 160 161 // convert boundary points to row/column coordinates 162 Matrix b (bound.size (), 2); 163 for (unsigned int i = 0; i < bound.size (); i++) 164 { 165 const int point = bound [i]; 166 b (i, ROW) = point % rows; 167 b (i, COL) = point / rows; 168 } 169 170 retval.append (b); 171 return retval; 172 } 173