Polygon width and centerline

Hello everybody, i'm dealing with a head breaker idea...

i have a .shp polygon with a elongated shape like a snake, with variable width (like a snake that have eat an elephant).

I`m trying to obtain the centerlines and the width of the snake polygon at each point in the center line. I been thinking to slice the snake polygon into various sub polygon to extrac the with on each one and determine the centerline. But i'm lazzy.

Not i'm thinking that the better way could be extract a continuous data of the half width of the polygon along the centerline, not using the boring profiles, but using a matrix (raster).

Recently i have transformed the snake polygon to a snake raster in arcgis, so i have a raster area in .asc format with only a value (0) for the snake area and NoData for the rest.

I'm now trying to calculate buffer using this script:

http://www.mathworks.com/matlabcentral/fileexchange/11095-bufferm2

But i don't now how to continue, i must create a complement of the snake raster to calculate the buffers maybe?. But i'm stucked... i dont know how.

Somebody can bring me an idea or maybe a link to other guys with a similar issue?

Answer 1

[Mark Hayworth]

Easy. Starting from a binary image...

1. Use bwmorph() to get the skeleton.
2. Use bwdist() to get the distance from the centerline to the edge.
3. Then multiply them to get the centerline where every point on the centerline is the distance to the edge.

It's like 3 lines of code.

[Jules Ray]

that sounds pretti interesting.. but i'm not so familiarized with this kind of processing in matlab... can you include this lines of code to try it...

and thanks a lot for you help

Answer 2

[Image Analyst]

It's literally only 3 simple lines, doing exactly what he said with the functions he told you to use. Nonetheless, here's a full blown demo for you. It's more than 3 lines because I put in comments and stuff to read in the demo image and display the intermediate steps.

clc;    % Clear the command window.
close all;  % Close all figures (except those of imtool.)
imtool close all;  % Close all imtool figures.
clear;  % Erase all existing variables.
workspace;  % Make sure the workspace panel is showing.
fontSize = 20;

% Check that user has the Image Processing Toolbox installed.
hasIPT = license('test', 'image_toolbox');
if ~hasIPT
  % User does not have the toolbox installed.
  message = sprintf('Sorry, but you do not seem to have the Image Processing Toolbox.\nDo you want to try to continue anyway?');
  reply = questdlg(message, 'Toolbox missing', 'Yes', 'No', 'Yes');
  if strcmpi(reply, 'No')
    % User said No, so exit.
    return;
  end
end

% Read in a standard MATLAB gray scale demo image.
folder = fullfile(matlabroot, '\toolbox\images\imdemos');
baseFileName = 'circles.png';
% Get the full filename, with path prepended.
fullFileName = fullfile(folder, baseFileName);
% Check if file exists.
if ~exist(fullFileName, 'file')
  % File doesn't exist -- didn't find it there.  Check the search path for it.
  fullFileName = baseFileName; % No path this time.
  if ~exist(fullFileName, 'file')
    % Still didn't find it.  Alert user.
    errorMessage = sprintf('Error: %s does not exist in the search path folders.', fullFileName);
    uiwait(warndlg(errorMessage));
    return;
  end
end
binaryImage = imread(fullFileName);
% Fill Holes
binaryImage = imfill(binaryImage, 'holes');
% Display the original gray scale image.
subplot(2, 2, 1);
imshow(binaryImage, []);
title('Original Grayscale Image', 'FontSize', fontSize);
% Enlarge figure to full screen.
set(gcf, 'units','normalized','outerposition',[0 0 1 1]);
% Give a name to the title bar.
set(gcf,'name','Demo by ImageAnalyst','numbertitle','off') 

% Now we're ready to begin, since we have our binary image.
% Skeletonize it:
skeletonizedImage = bwmorph(binaryImage, 'skel', inf);
% Display it.
subplot(2, 2, 2);
imshow(skeletonizedImage, []);
title('Skeletonized Image', 'FontSize', fontSize);

% Now get the distance transform.
edtImage = bwdist(~binaryImage);
% Display it.
subplot(2, 2, 3);
imshow(edtImage, []);
title('Euclidean Distance Transform', 'FontSize', fontSize);

% Now multiply
distanceFromEdge = edtImage .* single(skeletonizedImage);
% Display it.
subplot(2, 2, 4);
imshow(distanceFromEdge, []);
title('Distance From Edge', 'FontSize', fontSize);

[Jules Ray]

Thnks i gonna try with this code.. looks good

but a question: only works with .png files or .jpg.... what about a raster image like a .asc (with 1 and 0 gridcodes)?

thanks for your anwer

[Image Analyst]

It will work with any image as long as you can get it into a 2D array variable. I've never heard of an asc file. If it's a text file you might have to use dlmread() or importdata() instead of imread().

[Jules Ray]

Thnk you very much... it`s working amazingly good... even with my .asc (ascci table)

just to general knowledge .asc is an ascii file, in fact a table...