identify edges of a boundary THIS IS NO USER FUNCTION
0001 function varargout = identifyEdges(bnd, optionStruct) 0002 % identify edges of a boundary 0003 % 0004 % THIS IS NO USER FUNCTION 0005 0006 % The elk-library: convex geometry applied to crystallization modeling. 0007 % Copyright (C) 2013 Alexander Reinhold 0008 % 0009 % This program is free software: you can redistribute it and/or modify it 0010 % under the terms of the GNU General Public License as published by the 0011 % Free Software Foundation, either version 3 of the License, or (at your 0012 % option) any later version. 0013 % 0014 % This program is distributed in the hope that it will be useful, but 0015 % WITHOUT ANY WARRANTY; without even the implied warranty of 0016 % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 0017 % General Public License for more details. 0018 % 0019 % You should have received a copy of the GNU General Public License along 0020 % with this program. If not, see <http://www.gnu.org/licenses/>. 0021 0022 % The main idea comes the Hough transform, an algorithm to identify lines 0023 % in a binary image. The idea is that every line has two parameters: 0024 % orientation and distance from origin. The Hough transfers an image with 0025 % pixely in x/y direction into an image with pixels according to 0026 % angles/distances of the edges. Hence, it discretizes the parameter 0027 % space. In principle, it takes each pixel and direction, determines the 0028 % required distance of the line and increases the corresponding bin of 0029 % the Hough-image. 0030 % In this approach, we will do the same for every sample point. Bright 0031 % points (maxima) in the Hough-image will appear for a line that fits the 0032 % boundary curve well, because lot's of points contribute to these 0033 % specific parameters. 0034 % To allow a reasonable detection of proper parameters, we also include the 0035 % understanding of roughnes and minimal edge lengths. If we assume a 0036 % certain allowed roughness, then much more points can be assigned to a 0037 % point in Hough-space. This increases the imporance of the maxima 0038 % compared to the lower values and allows to use the minimal edge length 0039 % in relation to the total number of pixels in the image. 0040 0041 %% Hough grid 0042 % basic grid spacing 0043 dMin = min(bnd.polarDistance)*0.95; 0044 dMax = max(bnd.polarDistance)*1.05; 0045 deltaAngle = optionStruct.tolAngle; 0046 deltaDistance = optionStruct.tolDistance; 0047 % grid points 0048 houghAngleVector = 0:deltaAngle:(2*pi); 0049 houghDistanceVector = dMin:deltaDistance:dMax; 0050 % grid point boundaries for distance 0051 distanceLowerBoundVector = houghDistanceVector - deltaDistance/2; 0052 distanceUpperBoundVector = houghDistanceVector + deltaDistance/2; 0053 % obtain grid (for later plotting) 0054 [dGrid phiGrid] = meshgrid(houghDistanceVector, houghAngleVector); 0055 % initialize hough image matrix 0056 houghMatrix = dGrid*0; 0057 0058 %% evaluate Hough values 0059 % normals according to angle grid points 0060 normalMatrix = computeNormalFromAngle(houghAngleVector); 0061 % determine distances according to angles 0062 distanceMatrix = normalMatrix' * [bnd.cartX; bnd.cartY]; 0063 % fill bins by cycling through available distances 0064 for iDist = 1:length(houghDistanceVector) 0065 thisPlus = sum(... 0066 (distanceMatrix < distanceUpperBoundVector(iDist)) & ... 0067 (distanceMatrix >= distanceLowerBoundVector(iDist)), 2); 0068 houghMatrix(:, iDist) = houghMatrix(:, iDist) + thisPlus; 0069 end 0070 % we standardize the values to the ratio of accumulated points by total 0071 % points 0072 houghMatrix = houghMatrix/length(bnd.cartX); 0073 0074 if optionStruct.debugHough > 1 0075 % hough plot 0076 figure 0077 contourf(phiGrid/2/pi*360, dGrid, ... 0078 houghMatrix * 100) 0079 colorbar 0080 title('accumulated counts in percent of total pixels') 0081 0082 % edge angle plot 0083 figure, hold on 0084 plot(houghAngleVector/2/pi*360, max(houghMatrix*100, [], 2)) 0085 title('maximum hough value in percent of total pixels') 0086 end 0087 0088 %% find maxima with angle-max from hough - sliding maxima 0089 [angleRatingVector distanceIndexVector] = max(houghMatrix, [], 2); 0090 edgeAngleIndexVector = identifyMaximaSlidingWindow(angleRatingVector, ... 0091 round(optionStruct.minDiffAngle / deltaAngle)); 0092 % remove obviously short edges 0093 % edgeAngleIndexVector(... 0094 % angleRatingVector(edgeAngleIndexVector) < ... 0095 % 0.3*max(angleRatingVector)) = []; 0096 0097 % assign H-representation 0098 % houghAngleVector(edgeAngleIndexVector)/pi 0099 hrep.A = computeNormalFromAngle(houghAngleVector(edgeAngleIndexVector))'; 0100 hrep.h = houghDistanceVector(distanceIndexVector(edgeAngleIndexVector))'; 0101 % tmp = distanceIndexVector(edgeAngleIndexVector) 0102 % houghDistanceVector(tmp) 0103 0104 % %% adjust facet distance to mean square 0105 % for iFacet = 1:size(hrep.A, 1) 0106 % distanceVector = hrep.A(iFacet, :) * [bnd.cartX; bnd.cartY] - ... 0107 % hrep.h(iFacet); 0108 % filter = abs(distanceVector) < 0.5*optionStruct.tolDistance; 0109 % hrep.h = hrep.h + mean(distanceVector(filter)); 0110 % end 0111 0112 varargout{1} = hrep; 0113 if nargout > 1 0114 hough.angleGrid = phiGrid; 0115 hough.distanceGrid = dGrid; 0116 hough.angleVector = houghAngleVector; 0117 hough.distanceVector = houghDistanceVector; 0118 hough.image = houghMatrix; 0119 varargout{2} = hough; 0120 end