summaryrefslogtreecommitdiffstats
path: root/SD-VBS/common/toolbox/toolbox_basic/affine/find_AD.m
diff options
context:
space:
mode:
authorLeo Chan <leochanj@live.unc.edu>2020-10-22 01:53:21 -0400
committerJoshua Bakita <jbakita@cs.unc.edu>2020-10-22 01:56:35 -0400
commitd17b33131c14864bd1eae275f49a3f148e21cf29 (patch)
tree0d8f77922e8d193cb0f6edab83018f057aad64a0 /SD-VBS/common/toolbox/toolbox_basic/affine/find_AD.m
parent601ed25a4c5b66cb75315832c15613a727db2c26 (diff)
Squashed commit of the sb-vbs branch.
Includes the SD-VBS benchmarks modified to: - Use libextra to loop as realtime jobs - Preallocate memory before starting their main computation - Accept input via stdin instead of via argc Does not include the SD-VBS matlab code. Fixes libextra execution in LITMUS^RT.
Diffstat (limited to 'SD-VBS/common/toolbox/toolbox_basic/affine/find_AD.m')
-rwxr-xr-xSD-VBS/common/toolbox/toolbox_basic/affine/find_AD.m82
1 files changed, 82 insertions, 0 deletions
diff --git a/SD-VBS/common/toolbox/toolbox_basic/affine/find_AD.m b/SD-VBS/common/toolbox/toolbox_basic/affine/find_AD.m
new file mode 100755
index 0000000..3cccefb
--- /dev/null
+++ b/SD-VBS/common/toolbox/toolbox_basic/affine/find_AD.m
@@ -0,0 +1,82 @@
1function [A,D] = find_AD(I,J,mask,w)
2%
3% [A,D] = find_AD(I,J,mask,w)
4%
5% find the matrix affine transform A and displacement D,
6% such that SSD difference of I(Ax-d)-J(x) is minimized,
7%
8
9%
10% Jianbo Shi
11%
12
13
14[gy1,gx1] = grad(I,w);
15[gy2,gx2] = grad(J,w);
16
17gx = 0.5*(gx1+gx2);
18gy = 0.5*(gy1+gy2);
19
20[size_y,size_x] = size(I);
21[center_y,center_x] = find_center(size_y,size_x);
22mask = mask(w+1:size_y-w,w+1:size_x-w);
23
24[x,y] = meshgrid(1:size_x,1:size_y);
25x = x- center_x;
26y = y-center_y;
27
28x = x(w+1:size_y-w,w+1:size_x-w);
29y = y(w+1:size_y-w,w+1:size_x-w);
30
31gx_sqr = gx.*mask.*gx;
32gx_gy = gx.*mask.*gy;
33gy_sqr = gy.*mask.*gy;
34
35x_sqr = x.*x;
36x_y = x.*y;
37y_sqr = y.*y;
38
39T= zeros(6,6);
40T(1,1) = 0.5*trapz(trapz(gx_sqr.*x_sqr));
41T(2,1) = trapz(trapz(gx_gy.*x_y));
42T(3,1) = trapz(trapz(gx_sqr.*x_y));
43T(4,1) = trapz(trapz(gx_gy.*x_sqr));
44T(5,1) = trapz(trapz(gx_sqr.*x));
45T(6,1) = trapz(trapz(gx_gy.*x));
46T(2,2) = 0.5*trapz(trapz(gy_sqr.*y_sqr));
47T(3,2) = trapz(trapz(gx_gy.*y_sqr));
48T(4,2) = trapz(trapz(gy_sqr.*x_y));
49T(5,2) = trapz(trapz(gx_gy.*y));
50T(6,2) = trapz(trapz(gy_sqr.*y));
51T(3,3) = 0.5*trapz(trapz(gx_sqr.*y_sqr));
52T(4,3) = trapz(trapz(gx_gy.*x_y));
53T(5,3) = trapz(trapz(gx_sqr.*y));
54T(6,3) = trapz(trapz(gx_gy.*y));
55T(4,4) = 0.5*trapz(trapz(gy_sqr.*x_sqr));
56T(5,4) = trapz(trapz(gx_gy.*x));
57T(6,4) = trapz(trapz(gy_sqr.*x));
58T(5,5) = 0.5*trapz(trapz(gx_sqr));
59T(6,5) = trapz(trapz(gx_gy));
60T(6,6) = 0.5*trapz(trapz(gy_sqr));
61
62T = T+T';
63
64J = J(w+1:size_y-w,w+1:size_x-w);
65I = I(w+1:size_y-w,w+1:size_x-w);
66
67
68diff = (J-I).*mask;
69b(1) = trapz(trapz(diff.*gx.*x));
70b(2) = trapz(trapz(diff.*gy.*y));
71b(3) = trapz(trapz(diff.*gx.*y));
72b(4) = trapz(trapz(diff.*gy.*x));
73b(5) = trapz(trapz(diff.*gx));
74b(6) = trapz(trapz(diff.*gy));
75
76a = inv(T)*b';
77
78A = [1+a(1), a(3);
79a(4),1+a(2)];
80
81D= [a(5),a(6)];
82