imfeature.m

有关matlab的电子书籍有一定的帮助希望有用

function outstats = imfeature(varargin)
%IMFEATURE Compute feature measurements for image regions.
%   STATS = IMFEATURE(L,MEASUREMENTS) computes a set of
%   measurements for each labeled region in the label matrix
%   L. Positive integer elements of L correspond to different
%   regions. For example, the set of elements of L equal to 1
%   corresponds to region 1; the set of elements of L equal to 2
%   corresponds to region 2; and so on. STATS is a structure
%   array of length max(L(:)). The fields of the structure array
%   denote different measurements for each region, as specified
%   by MEASUREMENTS.
%
%   MEASUREMENTS can be a comma-separated list of strings, a cell
%   array containing strings, the string 'all', or the string
%   'basic'. The set of valid measurement strings includes:
%
%     'Area'              'ConvexHull'    'EulerNumber'
%     'Centroid'          'ConvexImage'   'Extrema'
%     'BoundingBox'       'ConvexArea'    'EquivDiameter'
%     'MajorAxisLength'   'Image'         'Solidity'
%     'MinorAxisLength'   'FilledImage'   'Extent'
%     'Orientation'       'FilledArea'    'PixelList'
%     'Eccentricity'
%
%   Measurement strings are case insensitive and can be
%   abbreviated.
%
%   If MEASUREMENTS is the string 'all', then all of the above
%   measurements are computed. If MEASUREMENTS is not specified
%   or if it is the string 'basic', then these measurements are 
%   computed: 'Area', 'Centroid', and 'BoundingBox'.
%
%   STATS = IMFEATURE(L,MEASUREMENTS,N) specifies the type of
%   connectivity used in computing the 'FilledImage',
%   'FilledArea', and 'EulerNumber' measurements. N can have a
%   value of either 4 or 8, where 4 specifies 4-connected objects
%   and 8 specifies 8-connected objects; if the argument is
%   omitted, it defaults to 8.
%
%   Class Support
%   -------------
%   The input label matrix L can be of class double or uint8.
%
%   See also BWLABEL, ISMEMBER.

%   Steven L. Eddins, September 1997
%   Copyright 1993-1998 The MathWorks, Inc.  All Rights Reserved.
%   $Revision: 1.6 $  $Date: 1997/11/24 15:35:23 $

officialStats = {'Area'
                 'Centroid'
                 'BoundingBox'
                 'MajorAxisLength'
                 'MinorAxisLength'
                 'Eccentricity'
                 'Orientation'
                 'ConvexHull'
                 'ConvexImage'
                 'ConvexArea'
                 'Image'
                 'FilledImage'
                 'FilledArea'
                 'EulerNumber'
                 'Extrema'
                 'EquivDiameter'
                 'Solidity'
                 'Extent'
                 'PixelList'};
tempStats = {'PerimeterCornerPixelList'};
allStats = [officialStats ; tempStats];

[L, requestedStats, n, msg] = ParseInputs(officialStats, ...
        varargin{:});
if (~isempty(msg))
    error(msg);
end

if (isempty(L))
    numObjs = 0;
else
    numObjs = round(max(double(L(:))));
end

% Initialize the stats structure array.
numStats = length(allStats);
empties = cell(numStats, numObjs);
stats = cell2struct(empties, allStats, 1);

% Initialize the computedStats structure array.
zz = cell(numStats, 1);
for k = 1:numStats
    zz{k} = 0;
end
computedStats = cell2struct(zz, allStats, 1);

for k = 1:length(requestedStats)
    switch requestedStats{k}
        
    case 'Area'
        [stats, computedStats] = ComputeArea(L, stats, computedStats);
        
    case 'FilledImage'
        [stats, computedStats] = ComputeFilledImage(L,stats,computedStats,n);
        
    case 'FilledArea'
        [stats, computedStats] = ComputeFilledArea(L,stats,computedStats,n);
        
    case 'ConvexArea'
        [stats, computedStats] = ComputeConvexArea(L, stats, computedStats);
        
    case 'Centroid'
        [stats, computedStats] = ComputeCentroid(L, stats, computedStats);
        
    case 'EulerNumber'
        [stats, computedStats] = ComputeEulerNumber(L,stats,computedStats,n);
        
    case 'EquivDiameter'
        [stats, computedStats] = ComputeEquivDiameter(L, stats, computedStats);
        
    case 'Extrema'
        [stats, computedStats] = ComputeExtrema(L, stats, computedStats);
        
    case 'BoundingBox'
        [stats, computedStats] = ComputeBoundingBox(L, stats, computedStats);
        
    case {'MajorAxisLength', 'MinorAxisLength', 'Orientation', 'Eccentricity'}
        [stats, computedStats] = ComputeEllipseParams(L, stats, computedStats);
        
    case 'Solidity'
        [stats, computedStats] = ComputeSolidity(L, stats, computedStats);
        
    case 'Extent'
        [stats, computedStats] = ComputeExtent(L, stats, computedStats);
        
    case 'ConvexImage'
        [stats, computedStats] = ComputeConvexImage(L, stats, computedStats);
        
    case 'ConvexHull'
        [stats, computedStats] = ComputeConvexHull(L, stats, computedStats);
        
    case 'Image'
        [stats, computedStats] = ComputeImage(L, stats, computedStats);

    case 'PixelList'
        [stats, computedStats] = ComputePixelList(L, stats, computedStats);
        
    end
end

% Initialize the output stats structure array.
numStats = length(requestedStats);
empties = cell(numStats, numObjs);
outstats = cell2struct(empties, requestedStats, 1);

% Initialize the subsref structure.
s(1).type = '()';
s(1).subs = {};
s(2).type = '.';
s(2).subs = '';

% Copy only the requested stats into the output.
for k = 1:numObjs
    for p = 1:length(requestedStats)
        s(1).subs = {k};
        s(2).subs = {requestedStats{p}};
        
        % In normal MATLAB syntax, the line below is the same as:
        %    outstats(k).fieldname = stats(k).fieldname
        % where fieldname is the string contained in
        % requestedStats{p}.  If you don't give subsasgn an
        % output argument it changes its first input argument
        % in-place.
        subsasgn(outstats, s, subsref(stats, s));
    end
end
    
%%%
%%% ComputeArea
%%%
function [stats, computedStats] = ComputeArea(L, stats, computedStats)
%   The area is defined to be the number of pixels belonging to
%   the region.

if ~computedStats.Area
    computedStats.Area = 1;

    [stats, computedStats] = ComputePixelList(L, stats, computedStats);

    for k = 1:length(stats)
        stats(k).Area = size(stats(k).PixelList, 1);
    end
end

%%%
%%% ComputeEquivDiameter
%%%
function [stats, computedStats] = ComputeEquivDiameter(L, stats, computedStats)
%   Computes the diameter of the circle that has the same area as
%   the region.
%   Ref: Russ, The Image Processing Handbook, 2nd ed, 1994, page
%   511.

if ~computedStats.EquivDiameter
    computedStats.EquivDiameter = 1;
    
    [stats, computedStats] = ComputeArea(L, stats, computedStats);

    factor = 2/sqrt(pi);
    for k = 1:length(stats)
        stats(k).EquivDiameter = factor * sqrt(stats(k).Area);
    end
end

%%%
%%% ComputeFilledImage
%%%
function [stats, computedStats] = ComputeFilledImage(L,stats,computedStats,n)
%   Uses bwfill with the 'holes' option to fill holes in the
%   region.  The last argument, n, specifies the foreground
%   connectivity and may be either 4 or 8.

if ~computedStats.FilledImage
    computedStats.FilledImage = 1;
    
    [stats, computedStats] = ComputeImage(L, stats, computedStats);
    
    for k = 1:length(stats)
        stats(k).FilledImage = bwfill(stats(k).Image,'holes',n);
    end
end

%%%
%%% ComputeConvexArea
%%%
function [stats, computedStats] = ComputeConvexArea(L, stats, computedStats)
%   Computes the number of "on" pixels in ConvexImage.

if ~computedStats.ConvexArea
    computedStats.ConvexArea = 1;
    
    [stats, computedStats] = ComputeConvexImage(L, stats, computedStats);
    
    for k = 1:length(stats)
        stats(k).ConvexArea = sum(stats(k).ConvexImage(:));
    end
end

%%%
%%% ComputeFilledArea
%%%
function [stats, computedStats] = ComputeFilledArea(L,stats,computedStats,n)
%   Computes the number of "on" pixels in FilledImage.

if ~computedStats.FilledArea
    computedStats.FilledArea = 1;
    
    [stats, computedStats] = ComputeFilledImage(L,stats,computedStats,n);

    for k = 1:length(stats)
        stats(k).FilledArea = sum(stats(k).FilledImage(:));
    end
end

%%%
%%% ComputeConvexImage
%%%
function [stats, computedStats] = ComputeConvexImage(L, stats, computedStats)
%   Uses ROIPOLY to fill in the convex hull.

if ~computedStats.ConvexImage
    computedStats.ConvexImage = 1;
    
    [stats, computedStats] = ComputeConvexHull(L, stats, computedStats);
    [stats, computedStats] = ComputeBoundingBox(L, stats, computedStats);
    
    for k = 1:length(stats)
        M = stats(k).BoundingBox(4);
        N = stats(k).BoundingBox(3);
        hull = stats(k).ConvexHull;
        if (isempty(hull))
            stats(k).ConvexImage = logical(uint8(zeros(M,N)));
        else
            firstRow = stats(k).BoundingBox(2) + 0.5;
            firstCol = stats(k).BoundingBox(1) + 0.5;
            r = hull(:,2) - firstRow + 1;
            c = hull(:,1) - firstCol + 1;
            stats(k).ConvexImage = roipoly(M, N, c, r);
        end
    end
end

%%%
%%% ComputeCentroid
%%%
function [stats, computedStats] = ComputeCentroid(L, stats, computedStats)
%   [mean(r) mean(c)]

if ~computedStats.Centroid
    computedStats.Centroid = 1;
    
    [stats, computedStats] = ComputePixelList(L, stats, computedStats);
    
    for k = 1:length(stats)
        if (isempty(stats(k).PixelList))
            % The reason for the empty special case is that the call
            % to mean below returns a scalar NaN and we need a 1-by-2
            % NaN.
            stats(k).Centroid = [NaN NaN];
        else
            stats(k).Centroid = mean(stats(k).PixelList,1);
        end
    end
end

%%%
%%% ComputeEulerNumber
%%%
function [stats, computedStats] = ComputeEulerNumber(L,stats,computedStats,n)
%   Calls BWEULER on 'Image'; last argument specifies foreground
%   connectivity and may be 4 or 8.

if ~computedStats.EulerNumber
    computedStats.EulerNumber = 1;
    
    [stats, computedStats] = ComputeImage(L, stats, computedStats);
    
    for k = 1:length(stats)
        stats(k).EulerNumber = bweuler(stats(k).Image,n);
    end
end

%%%
%%% ComputeExtrema
%%%
function [stats, computedStats] = ComputeExtrema(L, stats, computedStats)
%   A 8-by-2 array; each row contains the x and y spatial
%   coordinates for these extrema:  leftmost-top, rightmost-top,
%   topmost-right, bottommost-right, rightmost-bottom, leftmost-bottom,
%   bottommost-left, topmost-left. 
%   reference: Haralick and Shapiro, Computer and Robot Vision
%   vol I, Addison-Wesley 1992, pp. 62-64.

if ~computedStats.Extrema
    computedStats.Extrema = 1;
    
    [stats, computedStats] = ComputePixelList(L, stats, computedStats);
    
    for k = 1:length(stats)
        pixelList = stats(k).PixelList;
        if (isempty(pixelList))
            stats(k).Extrema = zeros(8,2) + 0.5;
        else
            r = pixelList(:,2);
            c = pixelList(:,1);
            
            minR = min(r);
            maxR = max(r);
            minC = min(c);
            maxC = max(c);
            
            minRSet = find(r==minR);
            maxRSet = find(r==maxR);
            minCSet = find(c==minC);
            maxCSet = find(c==maxC);

            % Points 1 and 2 are on the top row.
            r1 = minR;
            r2 = minR;
            % Find the minimum and maximum column coordinates for
            % top-row pixels.
            tmp = c(minRSet);
            c1 = min(tmp);
            c2 = max(tmp);
            
            % Points 3 and 4 are on the right column.
            % Find the minimum and maximum row coordinates for
            % right-column pixels.
            tmp = r(maxCSet);