Mathematica中文论坛-非官方
标题:
豆粑粑 matlab find circle hist
[打印本页]
作者:
meatball1982
时间:
2017-7-10 19:56
标题:
豆粑粑 matlab find circle hist
本帖最后由 meatball1982 于 2017-7-10 20:00 编辑
XP的媳妇儿的问题,想想挺好玩的。试试。
有一个什么图,里面是一堆堆的圆。
原来的思路是手动找圆,然后测量半径。然后手动统计分布。注意,是手动,手动。我的天啊!
matlab里有一个imfindcircles函数,老好老好了。
把图中的圆找出来,这之前,要把圆的边界收缩一下,更加容易实现找到圆的过程。
clear all
clc
clf
%% outline
% find circle 4 Xiaopeng
% it costs at least 500 RMB
%% main %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% read image -------------------------------------
% ima=imread('../file_imgs/img_cir_orig.jpg');
ima=imread('../file_imgs/p1-1.png');
se = strel('ball',10,3);
eroded = imerode(ima,se); % splot the circles
im=rgb2gray(eroded); % change to gray
im(im<=82)=0;
im(im>82)=500; % 2 values
[cen,rad,met]=imfindcircles(im,[9 14],'ObjectPolarity','bright',...
'Sensitivity',0.973,'EdgeThreshold',0.02); % find circles
[hi,bi]=hist(rad,40); % hist radius
hi = hi./sum(hi); % percents
save Mat_tm.mat
%% fit distribution %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
load Mat_tm.mat
% fit gaussian -----------------------
x_bin=bi;
x_plot= linspace(8,16,200);
x1=linspace(min(x_bin),max(x_bin),200);
y1 = hi;
fx=@(b,x)1/(b(1)*sqrt(2*pi*b(2))).*exp(-(x-b(3)).^2/(2*b(2)^2));
b0=[2,11,0.5];
b=b0;
for l=1:20
% Solve nonlinear curve-fitting problems in least-squares sense.
% just a guess.
b=lsqcurvefit(fx,b,x_bin,y1);
% Nonlinear regression,
b=nlinfit(x_bin,y1,fx,b);
end
b1=b;
y_hi1 = fx(b1,x_plot); % for plot gaussian distribution
% fit log gaussian --------------------------
fx_log =@(b,x) 1./x ./(b(1)*sqrt(2*pi*b(2))) .* exp(-(log(x)-b(3)).^2/(2*b(2)^2));
b0_log0 = [ 2 11 0.5 ];
b=b0;
for l=1:20
% Solve nonlinear curve-fitting problems in least-squares sense.
% just a guess.
b=lsqcurvefit(fx_log,b,x_bin,y1);
% Nonlinear regression,
b=nlinfit(x_bin,y1,fx_log,b);
end
b1_log=b;
y_hi1_log = fx_log(b1_log,x_plot); % for plot gaussian log distribution
save Mat_fit.mat
%% plot results %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
load Mat_fit.mat
sbp_width=0.85;
sbp_heig=0.8;
n_row = 2;
n_col = 3;
[out_pos]=fun_mm_subplot_pos(n_row,n_col,sbp_width,sbp_heig);
h=figure(1)
set(h, 'Position', [100, 100, 1000, 550]);
% subplot(2,3,1)
ax=axes('position',out_pos(1,:));
imshow(ima)
ax=axes('position',out_pos(2,:));
imshow(eroded)
ax=axes('position',out_pos(3,:));
imshow(im)
hBright = viscircles(cen, rad,'EdgeColor','r','linewidth',0.5);
ax=axes('position',out_pos(4,:));
imshow(ima)
hBright = viscircles(cen, rad+10,'EdgeColor','r','linewidth',0.5);
ax=axes('position',out_pos(5,:));
plot(rad,'ko','markerfacecolor','k','markersize',3)
axis([ 0 length(rad) 9 14])
ax=axes('position',out_pos(6,:));
barh(bi,hi,'c','edgecolor','none')
hold on
plot(y_hi1,x_plot,'b-','linewidth',2)
plot(y_hi1_log,x_plot,'m-','linewidth',2)
legend('hist','gau','log')
axis([0 1.1*max(hi) 9 14])
复制代码
res4subplot.png
(370.49 KB, 下载次数: 1434)
下载附件
2017-7-10 19:59 上传
欢迎光临 Mathematica中文论坛-非官方 (http://ilovemathematica.com/)
Powered by Discuz! X3.2