[37] | 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)); |
---|