1 | %
|
---|
2 | % %designed for the good of the world by Philip Torr based on ideas contained in
|
---|
3 | % copyright Philip Torr and Microsoft Corp 2002
|
---|
4 | %
|
---|
5 |
|
---|
6 |
|
---|
7 | % /*
|
---|
8 | % *****************************************
|
---|
9 | % */This gives a new method for estimating F that I developed (see citations below) from 7 points:
|
---|
10 |
|
---|
11 | %
|
---|
12 | % thus we can get a one parameter family of solutions (modulo scale) F = F1 + a F2
|
---|
13 | % exploiting det | F | = 0
|
---|
14 | % this gives us 1,2 or 3 distinct solutions for a
|
---|
15 |
|
---|
16 | %MORE DETAILS IN:
|
---|
17 | % @phdthesis{Torr:thesis,
|
---|
18 | % author="Torr, P. H. S.",
|
---|
19 | % title="Outlier Detection and Motion Segmentation",
|
---|
20 | % school=" Dept. of Engineering Science, University of Oxford",
|
---|
21 | % year=1995}
|
---|
22 | %
|
---|
23 | %
|
---|
24 | % @article{Torr97c,
|
---|
25 | % author="Torr, P. H. S. and Murray, D. W. ",
|
---|
26 | % title="The Development and Comparison of Robust Methods for Estimating the Fundamental Matrix",
|
---|
27 | % journal="IJCV",
|
---|
28 | % volume = 24,
|
---|
29 | % number = 3,
|
---|
30 | % pages = {271--300},
|
---|
31 | % year=1997
|
---|
32 | % }
|
---|
33 |
|
---|
34 | %A is a 7 x 9 design matrix;
|
---|
35 | %big_result is a 3x9 matrix of putative F's
|
---|
36 |
|
---|
37 | % returns 0,1,2 or 3 the number of fits gained.
|
---|
38 |
|
---|
39 |
|
---|
40 | %
|
---|
41 | function [no_F, big_result] = torr_F_constrained_fit(x1,y1,x2,y2,m3)
|
---|
42 |
|
---|
43 | A(:,1) = x1(:).* x2(:);
|
---|
44 | A(:,2) = y1(:).* x2(:);
|
---|
45 | A(:,3) = m3* x2(:);
|
---|
46 |
|
---|
47 | A(:,4) = x1(:).* y2(:);
|
---|
48 | A(:,5) = y1(:).* y2(:);
|
---|
49 | A(:,6) = m3* y2(:);
|
---|
50 |
|
---|
51 | A(:,7) = x1(:).* m3;
|
---|
52 | A(:,8) = y1(:).* m3;
|
---|
53 | A(:,9) = m3* m3;
|
---|
54 |
|
---|
55 | no_F = 0;
|
---|
56 | sizeA = size(A);
|
---|
57 | if ((sizeA(1) ~= 7) | (sizeA(2) ~= 9))
|
---|
58 | error('A is wrong shape');
|
---|
59 | return;
|
---|
60 | end
|
---|
61 | % A
|
---|
62 | % disp(sprintf('Rank of A = %d',rank(A)));
|
---|
63 | % disp('Null space of A');
|
---|
64 | % null(A)
|
---|
65 |
|
---|
66 |
|
---|
67 | [U,S,V] = svd(A);
|
---|
68 |
|
---|
69 | %the two solutions
|
---|
70 | r1 = V(:,length(V));
|
---|
71 | r2 = V(:,length(V)-1);
|
---|
72 | % disp('Torr')
|
---|
73 | % r1
|
---|
74 | % r2
|
---|
75 |
|
---|
76 | p = torr_calc_cubic_coefs(r1, r2);
|
---|
77 | a_roots = roots(p);
|
---|
78 | %disp('Roots')
|
---|
79 | %a_roots
|
---|
80 |
|
---|
81 | % % we now want to chose a solution so that the constrtaint
|
---|
82 | % % | aF_1 + (1-a)F_2 | = 0 is satisfies --- this gives us a cubic in a
|
---|
83 | %
|
---|
84 | % %coefficients of c_1 a^3 + c_2 a^2 + c_3 a + c_4
|
---|
85 | % cubic = zeros(4,1);
|
---|
86 | %
|
---|
87 | % cubic = torr_accumulate_coeff(cubic, r1, r2,1,5,9, 1);
|
---|
88 | % cubic = torr_accumulate_coeff(cubic, r1, r2,1,6,8, 0);
|
---|
89 | % cubic = torr_accumulate_coeff(cubic, r1, r2,2,4,9, 0);
|
---|
90 | % cubic = torr_accumulate_coeff(cubic, r1, r2,2,6,7, 1);
|
---|
91 | % cubic = torr_accumulate_coeff(cubic, r1, r2,3,4,8, 1);
|
---|
92 | % cubic = torr_accumulate_coeff(cubic, r1, r2,3,5,7, 0);
|
---|
93 | %
|
---|
94 | % % /* cubic(0) x^3 + cubic(1) x^2 + cubic(2) x + cubic(3) etc */
|
---|
95 | % a_roots = roots(cubic);
|
---|
96 |
|
---|
97 | for i_root = 1:3
|
---|
98 | if isreal(a_roots(i_root))
|
---|
99 | % disp(sprintf('Real root #%d = %0.8g',i_root,a_roots(i_root)))
|
---|
100 | no_F = no_F + 1;
|
---|
101 | big_result(no_F,:) = r1(:)' * a_roots(i_root) + r2(:)' * (1 - a_roots(i_root));
|
---|
102 |
|
---|
103 | % test that it is the correct solution.
|
---|
104 | % simply det(big_result) = 0
|
---|
105 | % disp(sprintf('det(aF1 + (1-a)F2) = %0.8g',det(reshape(big_result(no_F,:),3,3)')))
|
---|
106 |
|
---|
107 | end
|
---|
108 | end
|
---|
109 |
|
---|
110 |
|
---|
111 |
|
---|
112 |
|
---|
113 |
|
---|
114 | %accumulates the coefficients for the
|
---|
115 | % function cubic_out = torr_accumulate_coeff(cubic,f1, f2, a, b, c, plus)
|
---|
116 | % temp(1) = f1(a) * f1(b) * f1(c);
|
---|
117 | % temp(2) = f1(a) * f1(b) * f2(c);
|
---|
118 | % temp(3) = f1(a) * f2(b) * f1(c);
|
---|
119 | % temp(4) = f1(a) * f2(b) * f2(c);
|
---|
120 | % temp(5) = f2(a) * f1(b) * f1(c);
|
---|
121 | % temp(6) = f2(a) * f1(b) * f2(c);
|
---|
122 | % temp(7) = f2(a) * f2(b) * f1(c);
|
---|
123 | % temp(8) = f2(a) * f2(b) * f2(c);
|
---|
124 | %
|
---|
125 | % if plus
|
---|
126 | % cubic_out(1) = cubic(1) + temp(1) - temp(2) - temp(3) + temp(4) - temp(5) + temp(6) + temp(7) - temp(8);
|
---|
127 | % cubic_out(2) = cubic(2) + temp(2) + temp(3) - 2*temp(4) + temp(5) - 2*temp(6) - 2*temp(7)+ 3*temp(8);
|
---|
128 | % cubic_out(3) = cubic(3) +temp(4) + temp(6) + temp(7) - 3 * temp(8);
|
---|
129 | % cubic_out(4) = cubic(4) +temp(8);
|
---|
130 | % else
|
---|
131 | % cubic_out(1) = cubic(1) -temp(1) - temp(2) - temp(3) + temp(4) - temp(5) + temp(6) + temp(7) - temp(8);
|
---|
132 | % cubic_out(2) = cubic(2) -temp(2) + temp(3) - 2*temp(4) + temp(5) - 2*temp(6) - 2*temp(7)+ 3*temp(8);
|
---|
133 | % cubic_out(3) = cubic(3) -temp(4) + temp(6) + temp(7) - 3 * temp(8);
|
---|
134 | % cubic_out(4) = cubic(4) -temp(8);
|
---|
135 | % end
|
---|