Extract smooth surface from mask image and expand/shrink
[V, F, xx, B] = vb_mask_to_surf_expand(B, 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.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 = 1 : 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/6/14
Copyright (C) 2011, ATR All Rights Reserved.
License : New BSD License(see VBMEG_LICENSE.txt)0001 function [Vx, Fx, xx, B] = vb_mask_to_surf_expand(B, Para) 0002 % Extract smooth surface from mask image and expand/shrink 0003 % [V, F, xx, B] = vb_mask_to_surf_expand(B, 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.Nvertex = 5000 : # of vertex 0012 % Para.Nloop = 200 : iteration number for fit boundary 0013 % Para.Nlast = 10 : iteration number for smoothing by spring model 0014 % Para.tangent_rate = 0.3 : spring constant 0015 % Para.mask_ratio = 0.5 : mask image force constant 0016 % Para.mask_threshold = 0.3 : mask threthold 0017 % Para.vstep = 1 : voxel size of mask image [mm] 0018 % Para.Radius = [ 6 -6 ] : morphological smoothing [mm] 0019 % = [ 6 ] : morphological expansion 0020 % = [] : No morphological operation 0021 % --- 膨張モデル 0022 % Para.tangent_rate = 0.3 : バネ力係数 0023 % Para.mask_ratio = 0.5 : マスク強度力係数 ( 強度>閾値:外向き) 0024 % Para.mask_threshold = 0.3 : マスク強度閾値 ( 強度<閾値:内向き) 0025 % 0026 % 合力 = tangent_rate * (平滑化バネ力) 0027 % + mask_ratio * (外向き法線方向) * (マスク強度 - 閾値) 0028 % --- バネモデル平滑化 0029 % 0030 % 合力 = tangent_rate * (平滑化バネ力) 0031 % 0032 % 2006/10/12 M. Sato 0033 % 2006/11/11 M. Sato 0034 % M. Sato 2007/6/14 0035 % 0036 % Copyright (C) 2011, ATR All Rights Reserved. 0037 % License : New BSD License(see VBMEG_LICENSE.txt) 0038 0039 if ~exist('Para','var'), Para = []; end 0040 0041 % --- subsampling step for mask image [mm] 0042 % 3D 画像処理のための間引きステップサイズ 0043 if ~isfield(Para,'vstep'), 0044 step = 1; 0045 Para.vstep = step; 0046 else 0047 step = Para.vstep; 0048 end; 0049 0050 if ~isfield(Para,'Nloop'), Para.Nloop = 200; end 0051 if ~isfield(Para,'Nlast'), Para.Nlast = 20; end 0052 if ~isfield(Para,'Nvertex'), Para.Nvertex = 5000; end 0053 if ~isfield(Para,'tangent_rate'), Para.tangent_rate = 0.3; end 0054 if ~isfield(Para,'Radius'),Para.Radius = [6 -6];end 0055 if ~isfield(Para,'plot_mode'), Para.plot_mode = 0; end 0056 0057 Dim = size(B); 0058 0059 % Binarize mask image 0060 ix = find( B(:) > 0 ); 0061 B = zeros(Dim); 0062 B(ix) = 1; 0063 0064 % 0065 % --- Apply morphology operation 0066 % 0067 B = vb_morphology_operation(B, Para.Radius, step); 0068 0069 % 0070 % --- 3次元データの平滑化 0071 % 0072 tic 0073 fprintf('smoothing\n') 0074 B = smooth3(B,'gaussian',3); 0075 vb_ptime(toc); 0076 0077 % 0078 %------- マスクパターンの境界面に膨張 0079 % 0080 tic 0081 fprintf('--- Surface expansion\n') 0082 0083 % マスクパターンの境界座標 0084 V = vb_get_boundary(B); 0085 % Scale back to original image 0086 V = V*step - step/2; 0087 X0 = fix(mean(V/step)); % 中心点 0088 0089 [Vx, Fx, xx] = vb_make_inner_sphere(V, Para.Nvertex);% 内接球面作成 0090 0091 if Para.plot_mode > 1 0092 dmax = 10; 0093 vb_plot_slice_surf( B, Vx/step, Fx, X0(1), 'x',[1 1],'r-',dmax); 0094 vb_plot_slice_surf( B, Vx/step, Fx, X0(2), 'y',[1 1],'r-',dmax); 0095 vb_plot_slice_surf( B, Vx/step, Fx, X0(3), 'z',[1 1],'r-',dmax); 0096 end 0097 0098 [Vx, Fx, xx] = vb_surf_smooth_fit(Vx,Fx,xx, B, Para);% 境界面に膨張 0099 0100 vb_ptime(toc) 0101 0102 if Para.plot_mode > 1 0103 dmax = 10; 0104 vb_plot_slice_surf( B, Vx/step, Fx, X0(1), 'x',[1 1],'r-',dmax); 0105 vb_plot_slice_surf( B, Vx/step, Fx, X0(2), 'y',[1 1],'r-',dmax); 0106 vb_plot_slice_surf( B, Vx/step, Fx, X0(3), 'z',[1 1],'r-',dmax); 0107 end 0108 0109 % 0110 %------- バネモデル による平滑化 0111 % 0112 tic 0113 fprintf('--- Surface smoothing by spring model\n') 0114 0115 Para.Nloop = Para.Nlast; % バネモデル平滑化繰り返し回数 0116 [Vx, Fx, xx] = vb_surf_smooth(Vx,Fx,xx,Para); 0117 0118 vb_ptime(toc) 0119 0120 fprintf('# of patches : %d\n',size(Fx,1)); 0121 fprintf('# of vertices: %d\n',size(Vx,1)); 0122 0123 if Para.plot_mode > 1 0124 dmax = 10; 0125 vb_plot_slice_surf( B, Vx/step, Fx, X0(1), 'x',[1 1],'r-',dmax); 0126 vb_plot_slice_surf( B, Vx/step, Fx, X0(2), 'y',[1 1],'r-',dmax); 0127 vb_plot_slice_surf( B, Vx/step, Fx, X0(3), 'z',[1 1],'r-',dmax); 0128 end 0129