豆粑粑 matlab find circle hist

meatball1982

XP的媳妇儿的问题，想想挺好玩的。试试。
有一个什么图，里面是一堆堆的圆。
原来的思路是手动找圆，然后测量半径。然后手动统计分布。注意，是手动，手动。我的天啊！
matlab里有一个imfindcircles函数，老好老好了。
把图中的圆找出来，这之前，要把圆的边界收缩一下，更加容易实现找到圆的过程。

1. clear all
   
2. clc
   
3. clf
   
4. 
   
5. %% outline
   
6. % find circle 4 Xiaopeng
   
7. % it costs at least 500 RMB
   
8. 
   
9. %% main %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   
10. % read image -------------------------------------
    
11. % ima=imread('../file_imgs/img_cir_orig.jpg');
    
12. ima=imread('../file_imgs/p1-1.png');
    
13. se = strel('ball',10,3);            
    
14. eroded = imerode(ima,se);            % splot the circles
    
15. im=rgb2gray(eroded);                 % change to gray
    
16. im(im<=82)=0;
    
17. im(im>82)=500;                       % 2 values
    
18. 
    
19. 
    
20. [cen,rad,met]=imfindcircles(im,[9 14],'ObjectPolarity','bright',...
    
21.     'Sensitivity',0.973,'EdgeThreshold',0.02); % find circles
    
22. 
    
23. [hi,bi]=hist(rad,40);                          % hist radius
    
24. hi = hi./sum(hi);                              % percents
    
25. 
    
26. save Mat_tm.mat
    
27. 
    
28. 
    
29. %% fit distribution %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    
30. load Mat_tm.mat
    
31. % fit gaussian      -----------------------
    
32. x_bin=bi;
    
33. x_plot= linspace(8,16,200);
    
34. x1=linspace(min(x_bin),max(x_bin),200);
    
35. y1 = hi;
    
36. 
    
37. fx=@(b,x)1/(b(1)*sqrt(2*pi*b(2))).*exp(-(x-b(3)).^2/(2*b(2)^2));
    
38. b0=[2,11,0.5];
    
39. b=b0;
    
40. for l=1:20
    
41.     % Solve nonlinear curve-fitting problems in least-squares sense.
    
42.     % just a guess.
    
43.     b=lsqcurvefit(fx,b,x_bin,y1);
    
44.     % Nonlinear regression,
    
45.     b=nlinfit(x_bin,y1,fx,b);
    
46.    
    
47. end
    
48. b1=b;           
    
49. y_hi1 = fx(b1,x_plot);          % for plot gaussian distribution
    
50. 
    
51. % fit log gaussian --------------------------
    
52. fx_log =@(b,x) 1./x ./(b(1)*sqrt(2*pi*b(2))) .* exp(-(log(x)-b(3)).^2/(2*b(2)^2));
    
53. b0_log0 = [ 2 11 0.5  ];
    
54. 
    
55. b=b0;
    
56. for l=1:20
    
57.     % Solve nonlinear curve-fitting problems in least-squares sense.
    
58.     % just a guess.
    
59.     b=lsqcurvefit(fx_log,b,x_bin,y1);
    
60.     % Nonlinear regression,
    
61.     b=nlinfit(x_bin,y1,fx_log,b);
    
62. end
    
63. b1_log=b;
    
64. 
    
65. y_hi1_log = fx_log(b1_log,x_plot); % for plot gaussian log distribution
    
66. 
    
67. save Mat_fit.mat
    
68. 
    
69. 
    
70. %% plot results %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    
71. load Mat_fit.mat
    
72. sbp_width=0.85;
    
73. sbp_heig=0.8;
    
74. n_row = 2;
    
75. n_col = 3;
    
76. [out_pos]=fun_mm_subplot_pos(n_row,n_col,sbp_width,sbp_heig);
    
77. h=figure(1)
    
78. set(h, 'Position', [100, 100, 1000, 550]);
    
79. % subplot(2,3,1)
    
80. ax=axes('position',out_pos(1,:));
    
81. imshow(ima)
    
82. ax=axes('position',out_pos(2,:));
    
83. imshow(eroded)
    
84. ax=axes('position',out_pos(3,:));
    
85. imshow(im)
    
86. hBright = viscircles(cen, rad,'EdgeColor','r','linewidth',0.5);
    
87. ax=axes('position',out_pos(4,:));
    
88. imshow(ima)
    
89. hBright = viscircles(cen, rad+10,'EdgeColor','r','linewidth',0.5);
    
90. ax=axes('position',out_pos(5,:));
    
91. plot(rad,'ko','markerfacecolor','k','markersize',3)
    
92. axis([ 0 length(rad) 9 14])
    
93. ax=axes('position',out_pos(6,:));
    
94. barh(bi,hi,'c','edgecolor','none')
    
95. hold on
    
96. plot(y_hi1,x_plot,'b-','linewidth',2)
    
97. plot(y_hi1_log,x_plot,'m-','linewidth',2)
    
98. legend('hist','gau','log')
    
99. axis([0 1.1*max(hi) 9 14])

复制代码