1 | %create_nullspace Create null-space for scene with perspective camera. |
---|
2 | % |
---|
3 | % Parameters: |
---|
4 | % depths .. matrix of zeros and ones meaning whether corresponding element |
---|
5 | % in M is scaled (i.e. multiplied by known perspective depth) |
---|
6 | |
---|
7 | function [nullspace, result] = create_nullspace(M, depths, central, opt) |
---|
8 | |
---|
9 | if nargin < 4, opt.trial_coef = 1; |
---|
10 | opt.threshold = .01; end |
---|
11 | |
---|
12 | I = ~isnan(M(1:3:end,:)); [m n] = size(I); show_mod = 10; use_maxtuples = 0; |
---|
13 | if opt.verbose, fprintf(1, 'Used 4-tuples (.=%d): ', show_mod); tic; end |
---|
14 | |
---|
15 | if central, cols_scaled(1:n) = 0; cols_scaled(find(I(central,:) > 0)) = 1; |
---|
16 | else, cols_scaled = []; end |
---|
17 | |
---|
18 | num_trials = round(opt.trial_coef * n); |
---|
19 | |
---|
20 | if opt.verbose, fprintf(1,'(Allocating memory...'); end |
---|
21 | nullspace(size(M,1),num_trials) = 0; % Memory allocation: at least this |
---|
22 | % because at least one column is |
---|
23 | % added per one 4-tuple. |
---|
24 | width = 0; |
---|
25 | if opt.verbose, fprintf(1, ')'); end |
---|
26 | |
---|
27 | tenth = .1; % because of the first % |
---|
28 | result.used = 0; result.failed = 0; |
---|
29 | for i = 1:num_trials |
---|
30 | % choose a 4/max-tuple |
---|
31 | cols = 1:n; |
---|
32 | rows = 1:m; |
---|
33 | |
---|
34 | cols_chosen = []; t=1; failed = 0; |
---|
35 | if central, |
---|
36 | scaled_ensured = 0; |
---|
37 | else |
---|
38 | scaled_ensured = 1; % trial version: no scale controling when cutting |
---|
39 | end |
---|
40 | for t = 1:4 |
---|
41 | % choose one column, cut useless parts etc. |
---|
42 | [c, cols] = random_element(cols); |
---|
43 | cols_chosen = [cols_chosen c]; |
---|
44 | |
---|
45 | % check just added column |
---|
46 | rows = intersect(rows, find(I(:,c) > 0)); |
---|
47 | |
---|
48 | if t < 4, |
---|
49 | [rows, cols, scaled_ensured] = cut_useless(I, cols_scaled, ... |
---|
50 | cols_chosen, rows, cols, 4-t, scaled_ensured); |
---|
51 | end |
---|
52 | |
---|
53 | if isempty(rows), failed = 1; break; end |
---|
54 | end |
---|
55 | |
---|
56 | if ~failed |
---|
57 | % use the 4/max-tuple |
---|
58 | d = depths(rows,cols_chosen); |
---|
59 | % see ``debug code'' in the comment lower |
---|
60 | |
---|
61 | rowsbig = k2i(rows); |
---|
62 | submatrix=[]; for j=1:length(cols_chosen) % 4, |
---|
63 | submatrix=[ submatrix ... |
---|
64 | spread_depths_col(M(rowsbig,cols_chosen(j)), d(:,j)) ]; end |
---|
65 | debug=1; if debug, if size(submatrix, 1)<=size(submatrix,2) & opt.verbose |
---|
66 | fprintf(1,'-'); end;end |
---|
67 | subnull = nulleps(submatrix,opt.threshold); %svd(submatrix) |
---|
68 | if size(subnull,2)>0 & ( use_maxtuples | ... |
---|
69 | size(submatrix,1) == size(submatrix,2) + size(subnull,2)) |
---|
70 | nulltemp = zeros(size(M,1),size(subnull,2)); |
---|
71 | nulltemp(rowsbig,:) = subnull; % * (length(rows)/m); % weighting |
---|
72 | if width+size(nulltemp,2) > size(nullspace,2) % Memory allocation: |
---|
73 | if opt.verbose, fprintf(1,'(Allocating memory...)'); end |
---|
74 | mean_added = width/i; |
---|
75 | nullspace(size(M,1), size(nullspace,2) ... |
---|
76 | + round(mean_added * (num_trials-i))) = 0; |
---|
77 | end |
---|
78 | nullspace(:, width+1 : width+size(nulltemp,2)) = nulltemp; |
---|
79 | width = width +size(nulltemp,2); |
---|
80 | result.used = result.used +1; |
---|
81 | if mod(result.used, show_mod)==0 & opt.verbose, fprintf(1,'.'); end |
---|
82 | end |
---|
83 | else |
---|
84 | result.failed = result.failed +1; |
---|
85 | end |
---|
86 | |
---|
87 | if i/num_trials > .1*tenth |
---|
88 | if opt.verbose, fprintf(1,'%d%%', tenth*10); end |
---|
89 | if tenth < 1, tenth=0; end |
---|
90 | tenth = tenth + 1; |
---|
91 | end |
---|
92 | end |
---|
93 | |
---|
94 | if size(nullspace,2) > width, % cut off unused allocated memory |
---|
95 | if opt.verbose, fprintf(1,'(Cutting unused memory...'); end |
---|
96 | nullspace = nullspace(:,1:width); |
---|
97 | if opt.verbose, fprintf(1, ')'); end |
---|
98 | end |
---|
99 | if opt.verbose, fprintf(1, ['(' num2str(toc) ' sec)\n']); end |
---|
100 | result.tried = i; |
---|
101 | |
---|
102 | |
---|
103 | %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
104 | % debug code |
---|
105 | % if size(d, 1) == 1, bla; end |
---|
106 | % if sum(d(1,:)==0) > 2, blabla; end |
---|
107 | % if use_maxtuples | ... |
---|
108 | % size(d,1)*3 > size(d,2) + sum(d(:)==0)... % i.e. size(submatrix,1)> |
---|
109 | % ... % size(submatrix,2) |
---|
110 | % - sum(cols_scaled(cols_chosen)==0) % i.e. subtract 1 per each |
---|
111 | % % column full of zeros |
---|
112 | %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
113 | |
---|
114 | function [rows, cols, scaled_ensured] = cut_useless( ... %) |
---|
115 | I, cols_scaled, ... % this is always same |
---|
116 | cols_chosen, rows, cols, demanded, scaled_ensured) |
---|
117 | |
---|
118 | if ~scaled_ensured |
---|
119 | % check scaled columns |
---|
120 | if length(rows) == 2, demanded_scaled = 3; demanded_rows = 2; |
---|
121 | else demanded_scaled = 2; demanded_rows = 3; end |
---|
122 | cols_scaled_chosen = sum(cols_scaled(cols_chosen) > 0); |
---|
123 | |
---|
124 | % if no unscaled are allowed, they must be all cut |
---|
125 | if demanded == demanded_scaled - cols_scaled_chosen, |
---|
126 | cols = intersect(cols, find(cols_scaled > 0)); |
---|
127 | scaled_ensured = 1; |
---|
128 | end |
---|
129 | else demanded_rows = 2; end |
---|
130 | |
---|
131 | % check columns |
---|
132 | cols = cols(find(sum(I(rows,cols)) >= demanded_rows)); |
---|
133 | |
---|
134 | % check rows |
---|
135 | rows = rows(find(sum(I(rows,cols)') >= demanded)); |
---|
136 | |
---|
137 | function [element, rest] = random_element(set) |
---|
138 | % Take a random element out of a set. |
---|
139 | element = set(random_int(1, length(set))); |
---|
140 | rest = setdiff(set, element); |
---|
141 | |
---|
142 | function y = random_int(from,to) |
---|
143 | y = floor(from + (1 + to - from)*rand); |
---|
144 | |
---|
145 | function [N,s] = nulleps(M,tol) |
---|
146 | % Find the nullspace of M. This is the regular nullspace, augmented by |
---|
147 | % vectors that aren't really in the nullspace, but have low singular |
---|
148 | % values associated with them. tol is the threshold on these singular values. |
---|
149 | [u,s,v] = svd(M); |
---|
150 | sigsvs = sum(diag(s)>tol); |
---|
151 | N = u(:,sigsvs+1:size(u,2)); |
---|