1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
|
function dist = maskdist(vol)
%
% dist=maskdist(vol)
%
% return the distance in each voxel towards the nearest label boundaries
%
% author: Qianqian Fang (q.fang at neu.edu)
%
% input:
% vol: a 2D or 3D array
%
% output:
% dist: an integer array, storing the distance, in voxel unit, towards
% the nearest boundary between two distinct non-zero voxels, the
% zero voxels in the domain and space outside of the array
% are also treated as a unique non-zero label. If the goal is to
% get the minimum distance measured from the center of the voxel,
% one should use (dist-0.5).
%
% example:
%
% a=ones(60,60,60);
% a(:,:,1:10)=2;
% a(:,:,11:20)=3;
% im=maskdist(a);
% imagesc(squeeze(im(:,30,:)))
%
% -- this function is part of iso2mesh toolbox (http://iso2mesh.sf.net)
%
if (isempty(vol))
error('input vol can not be empty');
end
vals = unique(vol(:));
if (length(vals) > 256)
error('it appears that your input is a gray-scale image, you must convert it to binary or labels first');
end
newvol = ones(size(vol) + 2) * max(vals) + 1;
if (ndims(vol) == 2)
newvol(2:end - 1, 2:end - 1) = vol;
elseif (ndims(vol) == 3)
newvol(2:end - 1, 2:end - 1, 2:end - 1) = vol;
end
vals(end + 1) = newvol(1, 1, 1);
vals(vals == 0) = [];
newvol(newvol == 0) = newvol(1, 1, 1);
dist = ones(size(newvol)) * inf;
for i = 1:length(vals(:))
vv = (newvol == vals(i));
vdist = bwdist(vv);
vdist(vdist == 0) = inf;
dist = min(dist, vdist);
end
if (ndims(vol) == 2)
dist = dist(2:end - 1, 2:end - 1);
elseif (ndims(vol) == 3)
dist = dist(2:end - 1, 2:end - 1, 2:end - 1);
end
|
