Expand surface by morphology and smooth surface by spring model
[V_out, F_out, xx_out, B] = vb_surf_smooth_expand(V, F, xx, Para)
--- Input & Output
V(n, 1:3) : vertex
F(j, 1:3) : patch index
xx(n, 1:3) : normal vector
B(NX,NY,NZ) : mask image : inside of surface is filled out
--- Parameter
Para : parameter structure ( = default value )
Para.Dim : Dimension of mask image
Para.Nvertex = 5000 : # of vertex
Para.Nloop = 200 : iteration number for fit boundary
Para.Nlast = 10 : iteration number for smoothing by spring model
Para.tangent_rate = 0.3 : spring constant
Para.mask_ratio = 0.5 : mask image force constant
Para.mask_threshold = 0.3 : mask threthold
Para.vstep = 2 : voxel size of mask image [mm]
Para.Radius = [ 6 -6 ] : morphological smoothing [mm]
= [ 6 ] : morphological expansion
= [] : No morphological operation
--- 膨張モデル
Para.tangent_rate = 0.3 : バネ力係数
Para.mask_ratio = 0.5 : マスク強度力係数 ( 強度>閾値:外向き)
Para.mask_threshold = 0.3 : マスク強度閾値 ( 強度<閾値:内向き)
合力 = tangent_rate * (平滑化バネ力)
+ mask_ratio * (外向き法線方向) * (マスク強度 - 閾値)
--- バネモデル平滑化
合力 = tangent_rate * (平滑化バネ力)
2006/10/12 M. Sato
2006/11/11 M. Sato
M. Sato 2007-3-16
Copyright (C) 2011, ATR All Rights Reserved.
License : New BSD License(see VBMEG_LICENSE.txt)0001 function [Vx, Fx, xx, B] = vb_surf_smooth_expand(V, F, xx, Para) 0002 % Expand surface by morphology and smooth surface by spring model 0003 % [V_out, F_out, xx_out, B] = vb_surf_smooth_expand(V, F, xx, Para) 0004 % --- Input & Output 0005 % V(n, 1:3) : vertex 0006 % F(j, 1:3) : patch index 0007 % xx(n, 1:3) : normal vector 0008 % B(NX,NY,NZ) : mask image : inside of surface is filled out 0009 % --- Parameter 0010 % Para : parameter structure ( = default value ) 0011 % Para.Dim : Dimension of mask image 0012 % Para.Nvertex = 5000 : # of vertex 0013 % Para.Nloop = 200 : iteration number for fit boundary 0014 % Para.Nlast = 10 : iteration number for smoothing by spring model 0015 % Para.tangent_rate = 0.3 : spring constant 0016 % Para.mask_ratio = 0.5 : mask image force constant 0017 % Para.mask_threshold = 0.3 : mask threthold 0018 % Para.vstep = 2 : voxel size of mask image [mm] 0019 % Para.Radius = [ 6 -6 ] : morphological smoothing [mm] 0020 % = [ 6 ] : morphological expansion 0021 % = [] : No morphological operation 0022 % --- 膨張モデル 0023 % Para.tangent_rate = 0.3 : バネ力係数 0024 % Para.mask_ratio = 0.5 : マスク強度力係数 ( 強度>閾値:外向き) 0025 % Para.mask_threshold = 0.3 : マスク強度閾値 ( 強度<閾値:内向き) 0026 % 0027 % 合力 = tangent_rate * (平滑化バネ力) 0028 % + mask_ratio * (外向き法線方向) * (マスク強度 - 閾値) 0029 % --- バネモデル平滑化 0030 % 0031 % 合力 = tangent_rate * (平滑化バネ力) 0032 % 0033 % 2006/10/12 M. Sato 0034 % 2006/11/11 M. Sato 0035 % M. Sato 2007-3-16 0036 % 0037 % Copyright (C) 2011, ATR All Rights Reserved. 0038 % License : New BSD License(see VBMEG_LICENSE.txt) 0039 0040 if ~exist('Para','var'), Para = []; end 0041 0042 if ~isfield(Para,'Dim') | isempty(Para.Dim) 0043 error('Size of mask image is not specified'); 0044 else 0045 Dim = Para.Dim; 0046 end 0047 % --- subsampling step for mask image [mm] 0048 % 3D 画像処理のための間引きステップサイズ 0049 if ~isfield(Para,'vstep'), 0050 step = 2; 0051 Para.vstep = step; 0052 else 0053 step = Para.vstep; 0054 end; 0055 0056 if ~isfield(Para,'Nloop'), Para.Nloop = 200; end 0057 if ~isfield(Para,'Nlast'), Para.Nlast = 20; end 0058 if ~isfield(Para,'Nvertex'), Para.Nvertex = 5000; end 0059 if ~isfield(Para,'tangent_rate'), Para.tangent_rate = 0.3; end 0060 if ~isfield(Para,'Radius'),Para.Radius = [6 -6];end 0061 if ~isfield(Para,'plot_mode'), Para.plot_mode = 0; end 0062 0063 0064 %%% DEBUG 0065 Debug_mode = 0; 0066 0067 % 0068 % Closed surface check 0069 % 0070 disp('--- Closed surface check'); 0071 omega = vb_solid_angle_check(V,F); 0072 fprintf('Solid angle = 4pi x %f\n',omega); 0073 if abs(omega-1)>1e-10, 0074 disp('Surface is not closed.'); 0075 return; 0076 else 0077 disp('Surface is closed.'); 0078 end 0079 0080 % 0081 %-------- 表面をマスクに変換 0082 % 0083 fprintf('--- Make mask image from surface\n') 0084 dmin = step/2; % 三角形分割長さ 0085 0086 B = vb_surf_to_mask(V, F, step, dmin, Dim); 0087 0088 % 0089 %------- 内部点を塗りつぶす 0090 % 0091 fprintf('--- Fill out inside the surface \n') 0092 0093 filval = 1; % 内部点マスクの値 0094 level = 0.5; % この値より小さい値を塗りつぶす 0095 0096 X0 = fix(mean(V/step)); % 中心点 0097 B = vb_flood_fill_3d(B, X0, filval, level); 0098 0099 if Debug_mode == 2, return; end 0100 % 0101 % --- Apply morphology erosion 0102 % 0103 B = vb_morphology_operation(B, Para.Radius, step); 0104 0105 if Debug_mode == 1, return; end 0106 % 0107 % --- 3次元データの平滑化 0108 % 0109 tic 0110 fprintf('smoothing\n') 0111 B = smooth3(B,'gaussian',3); 0112 vb_ptime(toc); 0113 0114 % 0115 %------- マスクパターンの境界面に膨張 0116 % 0117 tic 0118 fprintf('--- Surface expansion\n') 0119 0120 [Vx, Fx, xx] = vb_make_inner_sphere(V, Para.Nvertex);% 内接球面作成 0121 [Vx, Fx, xx] = vb_surf_smooth_fit(Vx,Fx,xx, B, Para);% 境界面に膨張 0122 0123 vb_ptime(toc) 0124 0125 if Para.plot_mode > 1 0126 dmax = 10; 0127 vb_plot_slice_surf( B, Vx/step, Fx, X0(1), 'x',[1 1],'r-',dmax); 0128 vb_plot_slice_surf( B, Vx/step, Fx, X0(2), 'y',[1 1],'r-',dmax); 0129 vb_plot_slice_surf( B, Vx/step, Fx, X0(3), 'z',[1 1],'r-',dmax); 0130 end 0131 0132 % 0133 %------- バネモデル による平滑化 0134 % 0135 tic 0136 fprintf('--- Surface smoothing by spring model\n') 0137 0138 Para.Nloop = Para.Nlast; % バネモデル平滑化繰り返し回数 0139 [Vx, Fx, xx] = vb_surf_smooth(Vx,Fx,xx,Para); 0140 0141 vb_ptime(toc) 0142 0143 fprintf('# of patches : %d\n',size(Fx,1)); 0144 fprintf('# of vertices: %d\n',size(Vx,1)); 0145 0146 if Para.plot_mode > 1 0147 dmax = 10; 0148 vb_plot_slice_surf( B, Vx/step, Fx, X0(1), 'x',[1 1],'r-',dmax); 0149 vb_plot_slice_surf( B, Vx/step, Fx, X0(2), 'y',[1 1],'r-',dmax); 0150 vb_plot_slice_surf( B, Vx/step, Fx, X0(3), 'z',[1 1],'r-',dmax); 0151 end 0152