|
沙发
楼主 |
发表于 2018-1-11 12:50:24
|
只看该作者
更新:
3D版本的
主程序
- clear all
- clc
- clf
- %% outline
- % gen rand circle, sphere in 3D
- % according : http://www.ilovematlab.cn/thread-526875-1-1.html
- % gen cen
- % gen rad with cen
- %
- %% main
- % % parameters
- % n =200;
- % b =10;
- % rmin=0.1;
- % rmax=0.3;
- %
- % % gen initial center
- % [ cen ] = fun_gen_cen_3d(n,b);
- %
- % % distance
- % D=pdist(cen);
- % D_mat = squareform(D);
- %
- %
- % % generate r ind for cal ri+rj distance
- % x1 = [1:n]; % ind
- % [ind1,ind2]=ndgrid(x1,x1);
- % out_ind=[ind2(:),ind1(:)];
- % ind_cho = out_ind(:,1)<out_ind(:,2);
- % out_ind(~ind_cho,:)=[];
- %
- % % gen initial r
- % r = fun_gen_r(n,D_mat,rmin,rmax);
- % R_mat=squareform(sum(r(out_ind)')');
- %
- % % judge
- % dis=(D_mat-R_mat)>0;
- % n_flg = n^2-n;
- %
- % n_cont = 1;
- % while sum(dis(:))~= n_flg | ... % dis
- % sum(cen(:,1)-r < 0 )> 0 | ... % left bound
- % sum(cen(:,2)-r < 0 )> 0 | ... % bottom bound
- % sum(cen(:,3)-r < 0 )> 0 | ... % z dir
- % sum(cen(:,1)+r > b )> 0 | ... % right bound
- % sum(cen(:,2)+r > b )> 0 | ... % up bound
- % sum(cen(:,3)+r > b )> 0
- %
- % n_cont=n_cont+1
- % [ cen ] = fun_gen_cen_3d(n,b); % new centers
- % r = fun_gen_r(n,D_mat,rmin,rmax); % new r
- % R_mat=squareform(sum(r(out_ind)')'); % ri+rj matrix
- % D = pdist(cen);
- % D_mat = squareform(D); % center matrix
- % dis=(D_mat-R_mat)>0;
- % end
- %
- %
- % save Mat_4plot.mat
- load Mat_4plot.mat
- %% output result
- hold on
- plot3(cen(:,1),cen(:,2),cen(:,3),'k.','markersize',3)
- axis([0 b 0 b 0 b])
- axis square
- ax = gca;
- grid on
- ax.XTick = [0:1:10];
- ax.YTick = [0:1:10];
- ax.ZTick = [0:1:10];
- % the=linspace(0,2*pi,100);
- % the=[the,0];
- %
- for k = 1: n
- [x_sp,y_sp,z_sp]=sphere(20);
- x_sp = r(k)*x_sp + cen(k,1);
- y_sp = r(k)*y_sp + cen(k,2);
- z_sp = r(k)*z_sp + cen(k,3);
- surf(x_sp,y_sp,z_sp,r(k)*ones(size(x_sp)),'edgecolor','none')
- % x_r = r(k)*x_sp+cen(k,1);
- % y_r = r(k)*sin(the)+cen(k,2);
- % plot(x_r,y_r,'-','linewidth',2)
- end
- box on
- axis equal
- colormap(jet)
- alpha(0.5)
- xlabel('x')
- ylabel('y')
- zlabel('z')
复制代码
3D的中心程序
- function [ cen ] = fun_gen_cen_3d(n,b)
- tm=randperm(b*b*b);
- ind = sort(tm(1:n));
- tm_0 = zeros(b,b,b);
- tm_0(tm(1:n))=1;
- [ind_x,ind_y,ind_z]=findND(tm_0);
- % cen_x = ceil(ind/b)+0.5-rand(1,n)*0.35;
- % cen_y = mod(ind,b) +0.5-rand(1,n)*0.35;
- cen_x = ind_x' -0.5+rand(1,n)*0.35;
- cen_y = ind_y' -0.5+rand(1,n)*0.35;
- cen_z = ind_z' -0.5+rand(1,n)*0.35;
- cen=[cen_x;
- cen_y;
- cen_z;]';
-
- cen = abs(cen);
复制代码
|
|