matlab - find out the orientation, length and radius of capped rectangular object -

i have image shown fig.1. trying fit binary image capped rectangular (fig.2) figure out:

1. the orientation (the angle between long axis , horizontal axis)
2. the length (l) , radius (r) of object. best way it? help.

my naive idea using least square fit find out these information found out there no equation capped rectangle. in matlab there function called rectangle can create capped rectangle seems plot purpose.

i solved 2 different ways , have notes on each approach below. each method varies in complexity need decide best trade application.

first approach: least-squares-optimization: here used unconstrained optimization through matlab's fminunc() function. take @ matlab's see options can set prior optimization. made simple choices approach working you.

in summary, setup model of capped rectangle function of parameters, l, w, , theta. can include r if wish don't think need that; examining continuity half-semi-circles @ each edge, think may sufficient let r = w, inspection of model geometry. reduces number of optimization parameters one.

i made model of capped rectangle using boolean layers, see cappedrectangle() function below. result, needed function calculate finite difference gradients of model respect l, w, , theta. if don't provide these gradients fminunc(), attempt estimate these found matlab's estimates didn't work application, provided own part of error function gets called fminunc() (see below).

i didn't have data right-clicked on image above , downloaded: 'arihm.png'

to read data did (creates variable cdata):

image = importdata('arihm.png'); vars = fieldnames(image); = 1:length(vars) assignin('base', vars{i}, image.(vars{i})); end 

then converted double type , "cleaned-up" data normalizing. note: pre-processing important optimization work properly, , may have been needed since didn't have raw data (as mentioned downloaded image webpage question):

data = im2double(cdata); data = data / max(data(:)); figure(1); imshow(data); % looks same image above 

now image sizes:

ny = size(data,1); nx = size(data,2); 

note #1: might consider adding center of capped rectangle, (xc,yc), optimization parameters. these degrees of freedom make difference in overall fitting results (see comment on final error function values below). didn't set here can follow approach used l, w, , theta, add functionality finite difference gradients. need setup capped rectangle model function of (xc,yc).

edit: out of curiosity added optimization on capped rectangle center, see results @ bottom.

note #2: "continuity" @ ends of capped rectangle, let r = w. if like, can later include r explicit optimization parameter following examples l, w, theta. might want have r1 , r2 @ each endpoint variables?

below arbitrary starting values used illustrate example optimization. i don't know how information have in application in general, should try provide best initial estimates can.

l = 25; w = l; theta = 90; params0 = [l w theta]; 

note different results based on initial estimates.

next display starting estimate (the cappedrectangle() function defined later):

caprect0 = reshape(cappedrectangle(params0,nx,ny),nx,ny); figure(2); imshow(caprect0); 

define anonymous function error metric (errorfunc() listed below):

f = @(x)errorfunc(x,data); % define several optimization parameters fminunc(): options = optimoptions(@fminunc,'gradobj','on','tolx',1e-3, 'display','iter'); % call optimizer: tic [x,fval,exitflag,output] = fminunc(f,params0,options); time = toc; disp(['convergence time (sec) = ',num2str(time)]); % results: disp(['l0 = ',num2str(l),'; ', 'l estimate = ', num2str(x(1))]); disp(['w0 = ',num2str(w),'; ', 'w estimate = ', num2str(x(2))]); disp(['theta0 = ',num2str(theta),'; ', 'theta estimate = ', num2str(x(3))]); caprectestimate = reshape(cappedrectangle(x,nx,ny),nx,ny); figure(3); imshow(caprectestimate); 

below output fminunc (for more details on each column see matlab's help):

 iteration f(x) step optimality cg-iterations 0 0.911579 0.00465 1 0.860624 10 0.00457 1 2 0.767783 20 0.00408 1 3 0.614608 40 0.00185 1 .... , on ... 15 0.532118 0.00488281 0.000962 0 16 0.532118 0.0012207 0.000962 0 17 0.532118 0.000305176 0.000962 0 

you can see final error metric values have not decreased relative starting value, indicates me model function doesn't have enough degrees of freedom "fit" data well, consider adding optimization parameters, e.g., image center, discussed earlier.

edit: added optimization on capped rectangle center, see results @ bottom.

now print results (using 2011 macbook pro):

convergence time (sec) = 16.1053 l0 = 25; l estimate = 58.5773 w0 = 25; w estimate = 104.0663 theta0 = 90; theta estimate = 36.9024 

and display results:

edit: exaggerated "thickness" of fitting results above because model trying fit data while keeping center fixed, resulting in larger values w. see updated results @ bottom.

you can see comparing data final estimate relatively simple model starts resemble data well.

you can go further , calculate error bars estimates setting own monte-carlo simulations check accuracy function of noise , other degrading factors (with known inputs can generate produce simulated data).

below model function used capped rectangle (note: way did image rotation kind of sketchy numerically , not robust finite-differences quick , dirty , gets going):

function result = cappedrectangle(params, nx, ny) [x,y] = meshgrid(-(nx-1)/2:(nx-1)/2,-(ny-1)/2:(ny-1)/2); l = params(1); w = params(2); theta = params(3); % units degrees r = w; % define r1 , r2 displaced rounded edges: x1 = x - l; x2 = x + l; r1 = sqrt(x1.^2+y.^2); r2 = sqrt(x2.^2+y.^2); % capped rectangle prior rotation (theta = 0): temp = double( (abs(x) <= l) & (abs(y) <= w) | (r1 <= r) | (r2 <= r) ); cappedrectanglerotated = im2double(imrotate(mat2gray(temp), theta, 'bilinear', 'crop')); result = cappedrectanglerotated(:); return 

and need error function called fminunc:

function [error, df_dx] = errorfunc(params,data) ny = size(data,1); nx = size(data,2); % anonymous function model: model = @(params)cappedrectangle(params,nx,ny); % least-squares error (analogous chi^2 in literature): f = @(x)sum( (data(:) - model(x) ).^2 ) / sum(data(:).^2); % scalar error: error = f(params); [df_dx] = finitediffgrad(f,params); return 

as function calculate finite difference gradients:

function [df_dx] = finitediffgrad(fun,x) n = length(x); x = reshape(x,n,1); % pick small delta, dx should experimented with: dx = norm(x(:))/10; % define array of dx values; h_array = dx*eye(n); df_dx = zeros(size(x)); f = @(x) feval(fun,x); % finite diff approx (use "centered difference" error o(h^2);) j = 1:n hj = h_array(j,:)'; df_dx(j) = ( f(x+hj) - f(x-hj) )/(2*dx); end return 

---

second approach: use regionprops()

as others have pointed out can use matlab's regionprops(). overall think work best tuning , checking insure doing expect. approach call (it lot simpler first approach!):

data = im2double(cdata); data = round(data / max(data(:))); s = regionprops(data, 'orientation', 'majoraxislength', ... 'minoraxislength', 'eccentricity', 'centroid'); 

and struct result s:

>> s s = centroid: [345.5309 389.6189] majoraxislength: 365.1276 minoraxislength: 174.0136 eccentricity: 0.8791 orientation: 30.9354 

this gives enough information feed model of capped rectangle. @ first glance seems way go, seems have mind set on approach (maybe first approach above).

anyway, below image of results (in red) overlaid on top of data can see looks quite good:

---

edit: couldn't myself, suspected including image center optimization parameter, better results obtained, went ahead , did check. sure enough, same starting estimates used earlier in least-squares estimation, here results:

iteration f(x) step optimality cg-iterations 0 0.911579 0.00465 1 0.859323 10 0.00471 2 2 0.742788 20 0.00502 2 3 0.530433 40 0.00541 2 ... , on ... 28 0.0858947 0.0195312 0.000279 0 29 0.0858947 0.0390625 0.000279 1 30 0.0858947 0.00976562 0.000279 0 31 0.0858947 0.00244141 0.000279 0 32 0.0858947 0.000610352 0.000279 0 

by comparison earlier values can see new least-square error values quite bit smaller when including image center, confirming suspected earlier (so no big surprise).

the updated estimates capped rectangle parameters thus:

convergence time (sec) = 96.0418 l0 = 25; l estimate = 89.0784 w0 = 25; w estimate = 80.4379 theta0 = 90; theta estimate = 31.614 

and relative image array center get:

xc = -22.9107 yc = 35.9257 

the optimization takes longer results improved seen visual inspection:

if performance issue may want consider writing own optimizer or first try tuning matlab's optimization parameters, perhaps using different algorithm options well; see optimization options above.

here code updated model:

function result = cappedrectangle(params, nx, ny) [x,y] = meshgrid(-(nx-1)/2:(nx-1)/2,-(ny-1)/2:(ny-1)/2); % extract params make code more readable: l = params(1); w = params(2); theta = params(3); % units degrees xc = params(4); % new param: image center in x yc = params(5); % new param: image center in y % shift coordinates image center: x = x-xc; y = y-yc; % define r = w constraint: r = w; % define r1 , r2 rounded edges: x1 = x - l; x2 = x + l; r1 = sqrt(x1.^2+y.^2); r2 = sqrt(x2.^2+y.^2); temp = double( (abs(x) <= l) & (abs(y) <= w) | (r1 <= r) | (r2 <= r) ); cappedrectanglerotated = im2double(imrotate(mat2gray(temp), theta, 'bilinear', 'crop')); result = cappedrectanglerotated(:); 

and prior calling fminunc() adjusted parameter list:

l = 25; w = l; theta = 90; % set image center 0 intial guess: xc = 0; yc = 0; params0 = [l w theta xc yc]; 

enjoy.