/
embed.m
99 lines (83 loc) · 2.21 KB
/
embed.m
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
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
function ff = embed(x, mask, varargin)
%|function ff = embed(x, mask, varargin)
%| embed x in nonzero elements of (logical) mask
%| in
%| x [np (L)] the "nonzero" pixels (lexicographically stacked)
%| mask [(Nd)] logical array, np = sum(mask)
%| option
%| '*dim' {0|1} 0: [(N) (L)] (default); 1: return [(N) *L]
%| out
%| ff [(Nd) (L)] viewable image(s)
%| if input is sparse, output is full double
%|
%| See also: masker
%|
%| Copyright 2000-9-16, Jeff Fessler, University of Michigan
if nargin == 1 && streq(x, 'test'), embed_test, return, end
if nargin < 2, help(mfilename), error(mfilename), end
if isempty(mask) % trick: handle empty mask quickly, ignoring '*dim'
ff = x;
return
end
arg.prod_dim = 0;
if nargin > 2
arg = vararg_pair(arg, varargin, 'subs', {'*dim', 'prod_dim'});
end
if ~islogical(mask), error 'mask must be logical', end
dimx = size(x);
cl = class(x);
if issparse(x)
cl = 'double';
end
pL = prod(dimx(2:end));
if islogical(x)
% cl = 'double';
ff = false([numel(mask) pL]); % [np *L]
else
ff = zeros([numel(mask) pL], cl); % [np *L]
end
%if is_pre_v7
% ff = zeros([numel(mask) pL]); % [np *L]
%else
% ff = zeros([numel(mask) pL], cl); % [np *L]
%end
if pL > 1
ff(mask(:),:) = reshapee(x, [], pL);
else
ff(mask(:),1) = x;
end
if ~arg.prod_dim
if ndims(mask) == 2 && size(mask,2) == 1 % 1d cases, 2008-12-14
% trick: omit '1' in 2nd dimension for 1d cases:
ff = reshape(ff, [size(mask,1) dimx(2:end)]); % [(Nd) (L)]
else
ff = reshape(ff, [size(mask) dimx(2:end)]); % [(Nd) (L)]
end
end
function embed_test
ig = image_geom('nx', 512, 'ny', 500, 'dx', 1);
ig.mask = ig.circ > 0;
ig.mask = conv2(double(ig.mask), ones(2), 'same') > 0;
x = [1:sum(ig.mask(:))]';
cpu etic
f1 = embed(x, ig.mask);
cpu etoc 'time for one'
f2 = embed([x x], ig.mask);
cpu etoc 'time for two'
jf_equal(f1, f2(:,:,2))
x = repmat(x, [1 2 3]);
f3 = embed(x, ig.mask);
jf_equal(f1, f3(:,:,2))
x = ig.unitv;
x = x(ig.mask);
f4 = ig.embed(x);
x = sparse(double(x));
% f5 = ig.embed(x); % no! trick in image_geom.m for sparse x
f5 = embed(x, ig.mask);
jf_equal(f4, f5)
% test 1d case
mask = logical([0 1 1 1 0]');
x = [1 2 3]';
jf_equal(embed(x,mask), [0 1 2 3 0]')
x = [1 2 3]' * [1 1]; % [np 2]
jf_equal(embed(x,mask), [0 1 2 3 0]'*[1 1])