1## Copyright (C) 1999 Andy Adler <adler@sce.carleton.ca> 2## Copyright (C) 2008 Søren Hauberg <soren@hauberg.org> 3## Copyright (C) 2015 Hartmut Gimpel <hg_code@gmx.de> 4## Copyright (C) 2015 Carnë Draug <carandraug@octave.org> 5## 6## This program is free software; you can redistribute it and/or 7## modify it under the terms of the GNU General Public License as 8## published by the Free Software Foundation; either version 3 of the 9## License, or (at your option) any later version. 10## 11## This program is distributed in the hope that it will be useful, but 12## WITHOUT ANY WARRANTY; without even the implied warranty of 13## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 14## General Public License for more details. 15## 16## You should have received a copy of the GNU General Public License 17## along with this program; if not, see 18## <http:##www.gnu.org/licenses/>. 19 20## -*- texinfo -*- 21## @deftypefn {Function File} {[@var{bw}, @var{thresh}] =} edge (@var{im}, @var{method}, @dots{}) 22## Find edges using various methods. 23## 24## The image @var{im} must be 2 dimensional and grayscale. The @var{method} 25## must be a string with the string name. The other input arguments are 26## dependent on @var{method}. 27## 28## @var{bw} is a binary image with the identified edges. @var{thresh} is 29## the threshold value used to identify those edges. Note that @var{thresh} 30## is used on a filtered image and not directly on @var{im}. 31## 32## @seealso{fspecial} 33## 34## @end deftypefn 35## @deftypefn {Function File} {} edge (@var{im}, @qcode{"Canny"}) 36## @deftypefnx {Function File} {} edge (@var{im}, @qcode{"Canny"}, @var{thresh}) 37## @deftypefnx {Function File} {} edge (@var{im}, @qcode{"Canny"}, @var{thresh}, @var{sigma}) 38## Find edges using the Canny method. 39## 40## @var{thresh} is two element vector for the hysteresis thresholding. 41## The lower and higher threshold values are the first and second elements 42## respectively. If it is a scalar value, the lower value is set to 43## @code{0.4 * @var{thresh}}. 44## 45## @var{sigma} is the standard deviation to be used on the Gaussian filter 46## that is used to smooth the input image prior to estimating gradients. 47## Defaults to @code{sqrt (2)}. 48## 49## @end deftypefn 50## @deftypefn {Function File} {} edge (@var{im}, @qcode{"Kirsch"}) 51## @deftypefnx {Function File} {} edge (@var{im}, @qcode{"Kirsch"}, @var{thresh}) 52## @deftypefnx {Function File} {} edge (@var{im}, @qcode{"Kirsch"}, @var{thresh}, @var{direction}) 53## @deftypefnx {Function File} {} edge (@var{im}, @qcode{"Kirsch"}, @var{thresh}, @var{direction}, @var{thinning}) 54## Find edges using the Kirsch approximation to the derivatives. 55## 56## Edge points are defined as points where the length of the gradient exceeds 57## a threshold @var{thresh}. 58## 59## @var{thresh} is the threshold used and defaults to twice the square root of the 60## mean of the gradient squared of @var{im}. 61## 62## @var{direction} is the direction of which the gradient is 63## approximated and can be @qcode{"vertical"}, @qcode{"horizontal"}, 64## or @qcode{"both"} (default). 65## 66## @var{thinning} can be the string @qcode{"thinning"} (default) or 67## @qcode{"nothinning"}. This controls if a simple thinning procedure 68## is applied to the edge image such that edge points also need to have 69## a larger gradient than their neighbours. The resulting "thinned" 70## edges are only one pixel wide. 71## 72## @end deftypefn 73## @deftypefn {Function File} {} edge (@var{im}, @qcode{"Lindeberg"}) 74## @deftypefnx {Function File} {} edge (@var{im}, @qcode{"Lindeberg"}, @var{sigma}) 75## Find edges using the the differential geometric single-scale edge 76## detector by Tony Lindeberg. 77## 78## @var{sigma} is the scale (spread of Gaussian filter) at which the edges 79## are computed. Defaults to @code{2}. 80## 81## This method does not use a threshold value and therefore does not return 82## one. 83## 84## @end deftypefn 85## @deftypefn {Function File} {} edge (@var{im}, @qcode{"LoG"}) 86## @deftypefnx {Function File} {} edge (@var{im}, @qcode{"LoG"}, @var{thresh}, @var{sigma}) 87## Find edges by convolving with the Laplacian of Gaussian (LoG) filter, 88## and finding zero crossings. 89## 90## Only zero crossings where the filter response is larger than @var{thresh} 91## are retained. @var{thresh} is automatically computed as 0.75*@math{M}, 92## where @math{M} is the mean of absolute value of LoG filter response. 93## 94## @var{sigma} sets the spread of the LoG filter. Default to @code{2}. 95## 96## @end deftypefn 97## @deftypefn {Function File} {} edge (@var{im}, @qcode{"Prewitt"}, @dots{}) 98## Find edges using the Prewitt approximation to the derivatives. 99## 100## This method is the same as Kirsch except a different approximation 101## gradient is used. 102## 103## @end deftypefn 104## @deftypefn {Function File} {} edge (@var{im}, @qcode{"Roberts"}) 105## @deftypefnx {Function File} {} edge (@var{im}, @qcode{"Roberts"}, @var{thresh}) 106## @deftypefnx {Function File} {} edge (@var{im}, @qcode{"Roberts"}, @var{thresh}, @var{thinning}) 107## Find edges using the Roberts approximation to the derivatives. 108## 109## This method is similar to Kirsch except a different approximation 110## gradient is used, and the default @var{thresh} is multiplied by sqrt(1.5). 111## In addition, there it does not accept the @var{direction} argument. 112## 113## @end deftypefn 114## @deftypefn {Function File} {} edge (@var{im}, @qcode{"Sobel"}, @dots{}) 115## Find edges using the Sobel approximation to the derivatives. 116## 117## This method is the same as Kirsch except a different approximation 118## gradient is used. 119## 120## @end deftypefn 121## @deftypefn {Function File} {} edge (@var{im}, @qcode{"zerocross"}, @var{thresh}, @var{filter}) 122## Find edges by convolving with a user-supplied filter and finding zero 123## crossings. 124## 125## @end deftypefn 126## @deftypefn {Function File} {} edge (@var{im}, @qcode{"Andy"}, @var{thresh}, @var{params}) 127## Find edges by the original (Andy Adlers) @code{edge} algorithm. 128## The steps are 129## @enumerate 130## @item 131## Do a Sobel edge detection and to generate an image at 132## a high and low threshold. 133## @item 134## Edge extend all edges in the LT image by several pixels, 135## in the vertical, horizontal, and 45 degree directions. 136## Combine these into edge extended (EE) image. 137## @item 138## Dilate the EE image by 1 step. 139## @item 140## Select all EE features that are connected to features in 141## the HT image. 142## @end enumerate 143## 144## The parameters for the method is given in a vector: 145## @table @asis 146## @item params(1)==0 or 4 or 8 147## Perform x connected dilatation (step 3). 148## @item params(2) 149## Dilatation coefficient (threshold) in step 3. 150## @item params(3) 151## Length of edge extention convolution (step 2). 152## @item params(4) 153## Coefficient of extention convolution in step 2. 154## @end table 155## defaults = [8, 1, 3, 3] 156## 157## @end deftypefn 158 159function [bw_out, thresh, varargout] = edge (im, method = "sobel", varargin) 160 161 if (nargin < 1 || nargin > 5) 162 print_usage (); 163 endif 164 165 if (! isgray (im)) 166 error ("edge: IM must be a grayscale image"); 167 endif 168 169 method = tolower (method); 170 switch (method) 171 case {"sobel", "prewitt", "kirsch"} 172 [bw, thresh] = edge_kirsch_prewitt_sobel (im, method, varargin{:}); 173 case "roberts" 174 [bw, thresh, varargout{1:2}] = edge_roberts (im, varargin{:}); 175 case "log" 176 [bw, thresh] = edge_log (im, varargin{:}); 177 case "zerocross" 178 [bw, thresh] = edge_zerocross (im, varargin{:}); 179 case "canny" 180 [bw, thresh] = edge_canny (im, varargin{:}); 181 case "lindeberg" 182 [bw] = edge_lindeberg (im, varargin{:}); 183 case "andy" 184 [bw, thresh] = edge_andy (im, varargin{:}); 185 otherwise 186 error ("edge: unsupported edge detector `%s'", method); 187 endswitch 188 189 if (nargout == 0) 190 imshow (bw); 191 else 192 bw_out = bw; 193 endif 194endfunction 195 196function [bw, thresh] = edge_kirsch_prewitt_sobel (im, method, varargin) 197 if (numel (varargin) > 3) 198 error ("edge: %s method takes at most 3 extra arguments", method); 199 endif 200 201 ## This supports thinning, direction, and thresh arguments in any 202 ## order. Matlab madness I tell you. 203 varargin = tolower (varargin); 204 direction_mask = (strcmp (varargin, "both") 205 | strcmp (varargin, "horizontal") 206 | strcmp (varargin, "vertical")); 207 thinning_mask = (strcmp (varargin, "thinning") 208 | strcmp (varargin, "nothinning")); 209 thresh_mask = cellfun (@(x) isnumeric (x) && (isscalar (x) || isempty (x)), 210 varargin); 211 212 if (! all (direction_mask | thinning_mask | thresh_mask)) 213 error ("edge: %s method takes only THRESH, DIRECTION, and THINNING arguments", 214 method); 215 endif 216 217 if (nnz (direction_mask) > 1) 218 error ("edge: more than 1 direction argument defined"); 219 elseif (any (direction_mask)) 220 direction = varargin{direction_mask}; 221 else 222 direction = "both"; 223 endif 224 225 if (nnz (thinning_mask) > 1) 226 error ("edge: more than 1 thinning argument defined"); 227 elseif (any (thinning_mask) 228 && strcmp (varargin{thinning_mask}, "nothinning")) 229 thinning = false; 230 else 231 thinning = true; 232 endif 233 234 if (nnz (thresh_mask) > 1) 235 error ("edge: more than 1 threshold argument defined"); 236 elseif (any (thresh_mask)) 237 thresh = varargin{thresh_mask}; 238 else 239 thresh = []; 240 endif 241 242 ## For better speed, the calculation is done with 243 ## squared edge strenght values (strength2) and 244 ## squared threshold values (thresh2). 245 246 h1 = fspecial (method); 247 h1 ./= sum (abs (h1(:))); # normalize h1 248 im = im2double (im); 249 switch (direction) 250 case "horizontal" 251 strength2 = imfilter (im, h1, "replicate").^2; 252 case "vertical" 253 strength2 = imfilter (im, h1', "replicate").^2; 254 case "both" 255 strength2 = (imfilter (im, h1, "replicate").^2 256 + imfilter (im, h1', "replicate").^2); 257 otherwise 258 error ("edge: unknown DIRECTION `%s' for %s method", direction, method); 259 endswitch 260 261 if (isempty (thresh)) 262 thresh2 = 4 * mean (strength2(:)); 263 thresh = sqrt (thresh2); 264 else 265 thresh2 = thresh .^ 2; 266 endif 267 268 if (thinning) 269 ## Keep edge strengths for use in non-maximum 270 ## suppresion in the simple_thinning step. 271 strength2(strength2<=thresh2) = 0; 272 bw = simple_thinning (strength2); 273 else 274 bw = strength2 > thresh2; 275 endif 276endfunction 277 278function [bw, thresh, g45, g135] = edge_roberts (im, varargin) 279 if (numel (varargin) > 2) 280 error ("edge: Roberts method takes at most 2 extra arguments"); 281 endif 282 283 ## This supports thinning, direction, and thresh arguments in any 284 ## order. Matlab madness I tell you. 285 varargin = tolower (varargin); 286 thinning_mask = (strcmp (varargin, "thinning") 287 | strcmp (varargin, "nothinning")); 288 thresh_mask = cellfun (@(x) isnumeric (x) && (isscalar (x) || isempty (x)), 289 varargin); 290 291 if (! all (thinning_mask | thresh_mask)) 292 error ("edge: roberts method takes only THRESH and THINNING arguments"); 293 endif 294 295 if (nnz (thinning_mask) > 1) 296 error ("edge: more than 1 thinning argument defined"); 297 elseif (any (thinning_mask) 298 && strcmp (varargin{thinning_mask}, "nothinning")) 299 thinning = false; 300 else 301 thinning = true; 302 endif 303 304 if (nnz (thresh_mask) > 1) 305 error ("edge: more than 1 threshold argument defined"); 306 elseif (any (thresh_mask)) 307 thresh = varargin{thresh_mask}; 308 else 309 thresh = []; 310 endif 311 312 h1 = [1 0; 0 -1] ./ 2; 313 h2 = [0 1; -1 0] ./ 2; 314 g45 = imfilter (im, h1, "replicate"); 315 g135 = imfilter (im, h2, "replicate"); 316 strength2 = g45.^2 + g135.^2; 317 318 if (isempty (thresh)) 319 thresh2 = 6 * mean (strength2(:)); 320 thresh = sqrt (thresh2); 321 else 322 thresh2 = thresh .^ 2; 323 endif 324 325 if (thinning) 326 ## Keep edge strengths for use in non-maximum 327 ## suppresion in the simple_thinning step. 328 strength2(strength2<=thresh2) = 0; 329 bw = simple_thinning (strength2); 330 else 331 bw = strength2 > thresh2; 332 endif 333endfunction 334 335function [bw, thresh] = edge_log (im, thresh = [], sigma = 2) 336 if (nargin > 1 && (! isnumeric (thresh) || all (numel (thresh) != [0 1]))) 337 error ("edge: THRESH for LoG method must be a numeric scalar or empty"); 338 endif 339 if (nargin > 2 && (! isnumeric (sigma) || ! isscalar (sigma))) 340 error ("edge: SIGMA for LoG method must be a numeric scalar"); 341 endif 342 343 f = fspecial ("log", (2 * ceil (3*sigma)) +1, sigma); 344 g = conv2 (im, f, "same"); 345 346 if (isempty (thresh)) 347 thresh = 0.75*mean(abs(g(:))); 348 endif 349 350 bw = (abs (g) > thresh) & zerocrossings (g); 351endfunction 352 353function [bw, thresh] = edge_zerocross (im, thresh, f) 354 ## because filter is a required argument, so is thresh 355 if (nargin != 3) 356 error ("edge: a FILTER and THRESH are required for the zerocross method"); 357 elseif (! isnumeric (thresh) || all (numel (thresh) != [0 1])) 358 error ("edge: THRESH for zerocross method must be a numeric scalar or empty"); 359 elseif (! isnumeric (f)) 360 error ("edge: FILTER for zerocross method must be numeric"); 361 endif 362 363 g = conv2 (im, f, "same"); 364 if (isempty (thresh)) 365 thresh = mean (abs (g(:))); 366 endif 367 bw = (abs (g) > thresh) & zerocrossings(g); 368endfunction 369 370function [bw, thresh] = edge_canny (im, thresh = [], sigma = sqrt (2)) 371 if (nargin > 1 && (! isnumeric (thresh) || all (numel (thresh) != [0 1 2]))) 372 error ("edge: THRESH for Canny method must have 0, 1, or 2 elements"); 373 endif 374 if (nargin > 2 && (! isnumeric (sigma) || ! isscalar (sigma))) 375 error ("edge: SIGMA for Canny method must be a numeric scalar"); 376 endif 377 378 ## Gaussian filtering to change the edge scale. 379 ## Treat each dimensions separately for performance. 380 gauss = fspecial ("gaussian", [1 (8*ceil(sigma))], sigma); 381 im = im2double (im); 382 J = imfilter (im, gauss, "replicate"); 383 J = imfilter (J, gauss', "replicate"); 384 385 ## edge detection with Prewitt filter (treat dimensions separately) 386 p = [1 0 -1]/2; 387 Jx = imfilter (J, p, "replicate"); 388 Jy = imfilter (J, p', "replicate"); 389 Es = sqrt (Jx.^2 + Jy.^2); 390 Es_max = max (Es(:)); 391 if (Es_max > 0) 392 Es ./= Es_max; 393 endif 394 Eo = pi - mod (atan2 (Jy, Jx) - pi, pi); 395 396 if (isempty (thresh)) 397 tmp = mean(abs(Es(:))); 398 thresh = [0.4*tmp, tmp]; 399 elseif (numel (thresh) == 1) 400 thresh = [0.4*thresh thresh]; 401 else 402 thresh = thresh(:).'; # always return a row vector 403 endif 404 bw = nonmax_suppress(Es, Eo, thresh(1), thresh(2)); 405endfunction 406 407function [bw] = edge_lindeberg (im, sigma = 2) 408 if (nargin > 1 && (! isnumeric (sigma) || ! isscalar (sigma))) 409 error ("edge: SIGMA for Lindeberg method must be a numeric scalar"); 410 endif 411 412 ## Filters for computing the derivatives 413 Px = [-1 0 1; -1 0 1; -1 0 1]; 414 Py = [1 1 1; 0 0 0; -1 -1 -1]; 415 Pxx = conv2 (Px, Px, "full"); 416 Pyy = conv2 (Py, Py, "full"); 417 Pxy = conv2 (Px, Py, "full"); 418 Pxxx = conv2 (Pxx, Px, "full"); 419 Pyyy = conv2 (Pyy, Py, "full"); 420 Pxxy = conv2 (Pxx, Py, "full"); 421 Pxyy = conv2 (Pyy, Px, "full"); 422 ## Change scale 423 L = imsmooth (double (im), "Gaussian", sigma); 424 ## Compute derivatives 425 Lx = conv2 (L, Px, "same"); 426 Ly = conv2 (L, Py, "same"); 427 Lxx = conv2 (L, Pxx, "same"); 428 Lyy = conv2 (L, Pyy, "same"); 429 Lxy = conv2 (L, Pxy, "same"); 430 Lxxx = conv2 (L, Pxxx, "same"); 431 Lyyy = conv2 (L, Pyyy, "same"); 432 Lxxy = conv2 (L, Pxxy, "same"); 433 Lxyy = conv2 (L, Pxyy, "same"); 434 ## Compute directional derivatives 435 Lvv = Lx.^2.*Lxx + 2.*Lx.*Ly.*Lxy + Ly.^2.*Lyy; 436 Lvvv = Lx.^3.*Lxxx + 3.*Lx.^2.*Ly.*Lxxy ... 437 + 3.*Lx.*Ly.^2.*Lxyy + 3.*Ly.^3.*Lyyy; 438 ## Perform edge detection 439 bw = zerocrossings (Lvv) & Lvvv < 0; 440endfunction 441 442## The 'andy' edge detector that was present in older versions of 'edge'. 443function [imout, thresh] = edge_andy (im, thresh, param2) 444 [n,m]= size(im); 445 xx= 2:m-1; 446 yy= 2:n-1; 447 448 filt= [1 2 1;0 0 0; -1 -2 -1]/8; tv= 2; 449 imo= conv2(im, rot90(filt), 'same').^2 + conv2(im, filt, 'same').^2; 450 if nargin<2 || thresh==[]; 451 thresh= sqrt( tv* mean(mean( imo(yy,xx) )) ); 452 end 453# sum( imo(:)>thresh ) / prod(size(imo)) 454 dilate= [1 1 1;1 1 1;1 1 1]; tt= 1; sz=3; dt=3; 455 if nargin>=3 456 # 0 or 4 or 8 connected dilation 457 if length(param2) > 0 458 if param2(1)==4 ; dilate= [0 1 0;1 1 1;0 1 0]; 459 elseif param2(1)==0 ; dilate= 1; 460 end 461 end 462 # dilation threshold 463 if length(param2) > 2; tt= param2(2); end 464 # edge extention length 465 if length(param2) > 2; sz= param2(3); end 466 # edge extention threshold 467 if length(param2) > 3; dt= param2(4); end 468 469 end 470 fobliq= [0 0 0 0 1;0 0 0 .5 .5;0 0 0 1 0;0 0 .5 .5 0;0 0 1 0 0; 471 0 .5 .5 0 0;0 1 0 0 0;.5 .5 0 0 0;1 0 0 0 0]; 472 fobliq= fobliq( 5-sz:5+sz, 3-ceil(sz/2):3+ceil(sz/2) ); 473 474 xpeak= imo(yy,xx-1) <= imo(yy,xx) & imo(yy,xx) > imo(yy,xx+1) ; 475 ypeak= imo(yy-1,xx) <= imo(yy,xx) & imo(yy,xx) > imo(yy+1,xx) ; 476 477 imht= ( imo >= thresh^2 * 2); # high threshold image 478 imht(yy,xx)= imht(yy,xx) & ( xpeak | ypeak ); 479 imht([1,n],:)=0; imht(:,[1,m])=0; 480 481% imlt= ( imo >= thresh^2 / 2); # low threshold image 482 imlt= ( imo >= thresh^2 / 1); # low threshold image 483 imlt(yy,xx)= imlt(yy,xx) & ( xpeak | ypeak ); 484 imlt([1,n],:)=0; imlt(:,[1,m])=0; 485 486# now we edge extend the low thresh image in 4 directions 487 488 imee= ( conv2( imlt, ones(2*sz+1,1) , 'same') > tt ) | ... 489 ( conv2( imlt, ones(1,2*sz+1) , 'same') > tt ) | ... 490 ( conv2( imlt, eye(2*sz+1) , 'same') > tt ) | ... 491 ( conv2( imlt, rot90(eye(2*sz+1)), 'same') > tt ) | ... 492 ( conv2( imlt, fobliq , 'same') > tt ) | ... 493 ( conv2( imlt, fobliq' , 'same') > tt ) | ... 494 ( conv2( imlt, rot90(fobliq) , 'same') > tt ) | ... 495 ( conv2( imlt, flipud(fobliq) , 'same') > tt ); 496# imee(yy,xx)= conv2(imee(yy,xx),ones(3),'same') & ( xpeak | ypeak ); 497 imee= conv2(imee,dilate,'same') > dt; # 498 499% ff= find( imht==1 ); 500% imout = bwselect( imee, rem(ff-1, n)+1, ceil(ff/n), 8); 501 imout = imee; 502endfunction 503 504## An auxiliary function that performs a very simple thinning. 505## Strength is an image containing the edge strength. 506## bw contains a 1 in (r,c) if 507## 1) strength(r,c) is greater than both neighbours in the 508## vertical direction, OR 509## 2) strength(r,c) is greater than both neighbours in the 510## horizontal direction. 511## Note the use of OR. 512function bw = simple_thinning(strength) 513 [r c] = size(strength); 514 x = ( strength > [ zeros(r,1) strength(:,1:end-1) ] & ... 515 strength > [ strength(:,2:end) zeros(r,1) ] ); 516 y = ( strength > [ zeros(1,c); strength(1:end-1,:) ] & ... 517 strength > [ strength(2:end,:); zeros(1,c) ] ); 518 bw = x | y; 519endfunction 520 521## Auxiliary function. Finds the zero crossings of the 522## 2-dimensional function f. (By Etienne Grossmann) 523function z = zerocrossings(f) 524 z0 = f<0; ## Negative 525 [R,C] = size(f); 526 z = zeros(R,C); 527 z(1:R-1,:) |= z0(2:R,:); ## Grow 528 z(2:R,:) |= z0(1:R-1,:); 529 z(:,1:C-1) |= z0(:,2:C); 530 z(:,2:C) |= z0(:,1:C-1); 531 532 z &= !z0; ## "Positive zero-crossings"? 533endfunction 534 535 536## Test the madness of arbitrary order of input for prewitt, kirsch, 537## and sobel methods. 538%!test 539%! im = [ 540%! 249 238 214 157 106 69 60 90 131 181 224 247 252 250 250 541%! 250 242 221 165 112 73 62 91 133 183 225 248 252 250 251 542%! 252 246 228 173 120 78 63 90 130 181 224 248 253 251 251 543%! 253 248 232 185 132 87 62 80 116 170 217 244 253 251 252 544%! 253 249 236 198 149 101 66 71 101 155 206 238 252 252 252 545%! 254 250 240 210 164 115 73 69 92 143 196 232 252 253 252 546%! 70 70 68 61 49 36 24 22 26 38 52 63 70 70 70 547%! 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 548%! 62 63 62 59 51 42 33 25 22 26 36 45 56 60 62 549%! 252 253 252 246 221 190 157 114 90 90 118 157 203 235 248 550%! 251 253 254 251 233 209 182 136 103 92 107 139 185 225 245 551%! 251 253 254 253 243 227 206 163 128 108 110 133 175 217 242 552%! 252 253 254 254 249 241 228 195 164 137 127 139 172 212 239 553%! ] / 255; 554%! 555%! methods = {"kirsch", "prewitt", "sobel"}; 556%! for m_i = 1:numel (methods) 557%! method = methods{m_i}; 558%! 559%! bw = edge (im, method, 0.2, "both", "thinning"); 560%! assert (edge (im, method, 0.2), bw) 561%! 562%! args = perms ({0.2, "both", "thinning"}); 563%! for i = 1:rows (args) 564%! assert (edge (im, method, args{i,:}), bw) 565%! endfor 566%! 567%! bw = edge (im, method, 0.2, "vertical", "nothinning"); 568%! args = perms ({0.2, "vertical", "nothinning"}); 569%! for i = 1:rows (args) 570%! assert (edge (im, method, args{i,:}), bw) 571%! endfor 572%! 573%! bw = edge (im, method, 0.2, "vertical", "thinning"); 574%! args = perms ({0.2, "vertical"}); 575%! for i = 1:rows (args) 576%! assert (edge (im, method, args{i,:}), bw) 577%! endfor 578%! 579%! bw = edge (im, method, 0.2, "both", "nothinning"); 580%! args = perms ({0.2, "nothinning"}); 581%! for i = 1:rows (args) 582%! assert (edge (im, method, args{i,:}), bw) 583%! endfor 584%! endfor 585 586%!error <more than 1 threshold argument> 587%! bw = edge (rand (10), "sobel", 0.2, 0.4) 588%!error <more than 1 thinning argument> 589%! bw = edge (rand (10), "sobel", "thinning", "nothinning") 590%!error <more than 1 direction argument> 591%! bw = edge (rand (10), "sobel", "both", "both") 592%!error <only THRESH, DIRECTION, and THINNING arguments> 593%! bw = edge (rand (10), "sobel", [0.2 0.7], "both", "thinning") 594 595%!error <more than 1 threshold argument> 596%! bw = edge (rand (10), "kirsch", 0.2, 0.4) 597%!error <more than 1 thinning argument> 598%! bw = edge (rand (10), "kirsch", "thinning", "nothinning") 599%!error <more than 1 direction argument> 600%! bw = edge (rand (10), "kirsch", "both", "both") 601%!error <only THRESH, DIRECTION, and THINNING arguments> 602%! bw = edge (rand (10), "kirsch", [0.2 0.7], "both", "thinning") 603 604%!error <more than 1 threshold argument> 605%! bw = edge (rand (10), "prewitt", 0.2, 0.4) 606%!error <more than 1 thinning argument> 607%! bw = edge (rand (10), "prewitt", "thinning", "nothinning") 608%!error <more than 1 direction argument> 609%! bw = edge (rand (10), "prewitt", "both", "both") 610%!error <only THRESH, DIRECTION, and THINNING arguments> 611%! bw = edge (rand (10), "prewitt", [0.2 0.7], "both", "thinning") 612 613%!test 614%! im = [ 615%! 249 238 214 157 106 69 60 90 131 181 224 247 252 250 250 616%! 250 242 221 165 112 73 62 91 133 183 225 248 252 250 251 617%! 252 246 228 173 120 78 63 90 130 181 224 248 253 251 251 618%! 253 248 232 185 132 87 62 80 116 170 217 244 253 251 252 619%! 253 249 236 198 149 101 66 71 101 155 206 238 252 252 252 620%! 254 250 240 210 164 115 73 69 92 143 196 232 252 253 252 621%! 70 70 68 61 49 36 24 22 26 38 52 63 70 70 70 622%! 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 623%! 62 63 62 59 51 42 33 25 22 26 36 45 56 60 62 624%! 252 253 252 246 221 190 157 114 90 90 118 157 203 235 248 625%! 251 253 254 251 233 209 182 136 103 92 107 139 185 225 245 626%! 251 253 254 253 243 227 206 163 128 108 110 133 175 217 242 627%! 252 253 254 254 249 241 228 195 164 137 127 139 172 212 239 628%! ] / 255; 629%! 630%! bw = edge (im, "roberts", .2, "thinning"); 631%! assert (edge (im, "roberts", 0.2), bw) 632%! assert (edge (im, "roberts", "thinning", 0.2), bw) 633%! 634%! bw = edge (im, "roberts", .2, "nothinning"); 635%! assert (edge (im, "roberts", "nothinning", 0.2), bw) 636 637%!error <more than 1 threshold argument> 638%! bw = edge (rand (10), "roberts", 0.2, 0.4) 639%!error <more than 1 thinning argument> 640%! bw = edge (rand (10), "roberts", "thinning", "nothinning") 641%!error <only THRESH and THINNING arguments> 642%! bw = edge (rand (10), "roberts", "both", "thinning") 643 644## test how canny threshold arguments are always a row vector 645%!test 646%! im = rand (10); 647%! [~, thresh] = edge (im, "canny"); 648%! assert (size (thresh), [1 2]) 649%! [~, thresh] = edge (im, "canny", [.2 .6]); 650%! assert (thresh, [.2 .6]) 651%! [~, thresh] = edge (im, "canny", [.2; .6]); 652%! assert (thresh, [.2 .6]) 653 654## test SOBEL edge detector 655%!test 656%! in = zeros (5); 657%! in(3,3) = 1; 658%! 659%! E = logical ([ 660%! 0 0 0 0 0 661%! 0 0 1 0 0 662%! 0 1 0 1 0 663%! 0 0 1 0 0 664%! 0 0 0 0 0]); 665%! assert (edge (in), E) 666%! assert (edge (uint8 (in.*100)), E) 667%! assert (edge (in, "sobel"), E) 668%! assert (edge (in, "sobel", 0), E) 669%! assert (edge (in, "sobel", 1), false (5)) 670%! 671%! [E, auto_thresh] = edge (in); 672%! assert (auto_thresh, 0.2449, 1e-4) 673%! 674%! V = logical([ 675%! 0 0 0 0 0 676%! 0 1 0 1 0 677%! 0 1 0 1 0 678%! 0 1 0 1 0 679%! 0 0 0 0 0]); 680%! assert (edge (in, "sobel", 0, "vertical"), V) 681%! 682%! H = logical ([ 683%! 0 0 0 0 0 684%! 0 1 1 1 0 685%! 0 0 0 0 0 686%! 0 1 1 1 0 687%! 0 0 0 0 0]); 688%! assert (edge (in, "sobel", 0, "horizontal"), H) 689%! 690%! V = false (5); 691%! V(3,2) = true; 692%! V(3,4) = true; 693%! assert (edge (in, "sobel", [], "vertical"), V) 694%! 695%! H = false (5); 696%! H(2,3) = true; 697%! H(4,3) = true; 698%! assert (edge (in, "sobel", [], "horizontal"), H) 699 700%!test 701%! A = ones (5); 702%! A(3, 3) = 0; 703%! expected = logical ([ 704%! 0 0 0 0 0 705%! 0 0 1 0 0 706%! 0 1 0 1 0 707%! 0 0 1 0 0 708%! 0 0 0 0 0]); 709%! assert (edge (A), expected) 710 711## test PREWITT edge detector 712%!test 713%! in = zeros (5); 714%! in(3, 3) = 1; 715%! 716%! E = logical ([ 717%! 0 0 0 0 0 718%! 0 1 0 1 0 719%! 0 0 0 0 0 720%! 0 1 0 1 0 721%! 0 0 0 0 0]); 722%! 723%! assert (edge (in, "prewitt"), E) 724%! 725%! [~, auto_thresh] = edge (in, "prewitt"); 726%! assert (auto_thresh, 0.2309, 1e-4) 727%! 728%! V = logical([ 729%! 0 0 0 0 0 730%! 0 1 0 1 0 731%! 0 1 0 1 0 732%! 0 1 0 1 0 733%! 0 0 0 0 0]); 734%! assert (edge (in, "prewitt", 0, "vertical"), V) 735%! 736%! H = logical ([ 737%! 0 0 0 0 0 738%! 0 1 1 1 0 739%! 0 0 0 0 0 740%! 0 1 1 1 0 741%! 0 0 0 0 0]); 742%! assert (edge (in, "prewitt", 0, "horizontal"), H) 743 744## test ROBERTS edge detector 745%!test 746%! in = zeros (5); 747%! in(3,3) = 1; 748%! in(3,4) = 0.9; 749%! 750%! E = logical ([ 751%! 0 0 0 0 0 752%! 0 0 1 0 0 753%! 0 0 1 0 0 754%! 0 0 0 0 0 755%! 0 0 0 0 0]); 756%! 757%! assert (edge (in, "roberts"), E) 758%! 759%! [~, auto_thresh] = edge (in, "roberts"); 760%! assert (auto_thresh, 0.6591, 1e-4) 761%! 762%! E45 = [0 0 0 0 0 763%! 0 -0.5 -0.45 0 0 764%! 0 0 0.50 0.45 0 765%! 0 0 0 0 0 766%! 0 0 0 0 0]; 767%! E135 = [0 0 0 0 0 768%! 0 0 -0.50 -0.45 0 769%! 0 0.5 0.45 0 0 770%! 0 0 0 0 0 771%! 0 0 0 0 0]; 772%! 773%! [~, ~, erg45, erg135] = edge (in, "roberts"); 774%! assert (erg45, E45) 775%! assert (erg135, E135) 776 777## test CANNY edge detector 778%!xtest 779%! ## The edge image is correct and Matlab compatible so those should 780%! ## pass. However, the threshold values used to generate the edge 781%! ## image are not the same as Matlab. 782%! 783%! in_8 = fspecial ("gaussian", [8 8], 2); 784%! in_8 /= in_8(4,4); 785%! in_8_uint8 = im2uint8 (in_8); 786%! 787%! ## Matlab changed their implementation of the Canny method in 788%! ## release 2011a. We are compatible with their new implementation 789%! ## but for testing purposes, this is the expected result for the 790%! ## old implementation. 791%! out_8_old = logical ([ 792%! 0 0 0 0 0 0 0 0 793%! 0 0 0 1 1 0 0 0 794%! 0 0 1 0 0 1 0 0 795%! 0 1 0 0 0 0 1 0 796%! 0 1 0 0 0 0 1 0 797%! 0 0 1 0 0 1 0 0 798%! 0 0 0 1 1 0 0 0 799%! 0 0 0 0 0 0 0 0]); 800%! 801%! out_8 = logical ([ 802%! 0 0 0 0 0 0 0 0 803%! 0 1 1 1 1 1 0 0 804%! 0 1 0 0 0 1 0 0 805%! 0 1 0 0 0 1 0 0 806%! 0 1 0 0 0 1 0 0 807%! 0 1 1 1 1 1 0 0 808%! 0 0 0 0 0 0 0 0 809%! 0 0 0 0 0 0 0 0]); 810%! out_thresh = [0.34375 0.859375]; 811%! 812%! [obs_edge, obs_thresh] = edge (in_8, "Canny"); 813%! assert (obs_edge, out_8) 814%! assert (obs_thresh, out_thresh) 815%! 816%! [obs_edge_givethresh, obs_thresh_givethresh] ... 817%! = edge (in_8, "Canny", out_thresh); 818%! assert (obs_edge_givethresh, out_8) 819%! assert (obs_thresh_givethresh, out_thresh) 820%! 821%! [obs_edge_uint8, obs_thresh_uint8] = edge (in_8_uint8, "Canny"); 822%! assert (obs_edge_uint8, out_8) 823%! assert (obs_thresh_uint8, out_thresh) 824 825%!xtest 826%! ## The edge image is correct and Matlab compatible so those should 827%! ## pass. However, the threshold values used to generate the edge 828%! ## image are not the same as Matlab. 829%! 830%! in_9 = fspecial ("gaussian", [9 9], 2); 831%! in_9 /= in_9(5,5); 832%! 833%! ## Matlab changed their implementation of the Canny method in 834%! ## release 2011a. We are compatible with their new implementation 835%! ## but for testing purposes, this is the expected result for the 836%! ## old implementation. 837%! out_9_old = logical ([ 838%! 0 0 0 0 0 0 0 0 0 839%! 0 0 0 0 0 0 0 0 0 840%! 0 0 0 1 1 1 0 0 0 841%! 0 0 1 0 0 0 1 0 0 842%! 0 0 1 0 0 0 1 0 0 843%! 0 0 1 0 0 0 1 0 0 844%! 0 0 0 1 1 1 0 0 0 845%! 0 0 0 0 0 0 0 0 0 846%! 0 0 0 0 0 0 0 0 0]); 847%! 848%! out_9 = logical ([ 849%! 0 0 0 0 0 0 0 0 0 850%! 0 0 1 1 1 1 0 0 0 851%! 0 1 1 0 0 1 1 0 0 852%! 0 1 0 0 0 0 1 0 0 853%! 0 1 0 0 0 0 1 0 0 854%! 0 1 1 0 0 1 1 0 0 855%! 0 0 1 1 1 1 0 0 0 856%! 0 0 0 0 0 0 0 0 0 857%! 0 0 0 0 0 0 0 0 0]); 858%! out_thresh = [0.35 0.875]; 859%! 860%! [obs_edge, obs_thresh] = edge (in_9, "Canny"); 861%! assert (obs_edge, out_9) 862%! assert (obs_thresh, out_thresh) 863%! 864%! [obs_edge_givethresh, obs_thresh_givethresh] ... 865%! = edge (in_9, "Canny", out_thresh); 866%! assert (obs_edge_givethresh, out_9) 867%! assert (obs_thresh_givethresh, out_thresh) 868