1%% Reset everything
2
3clear all;
4clc;
5close all;
6addpath('helpers');
7
8%% Configure the benchmark
9
10% noncentral case
11cam_number = 4;
12% Getting 10 points, and testing all algorithms with the respective number of points
13pt_number = 50;
14% noise test, so no outliers
15outlier_fraction = 0.0;
16% repeat 5000 tests per noise level
17iterations = 5000;
18
19% The algorithms we want to test
20algorithms = { 'seventeenpt'; 'rel_nonlin_noncentral'; 'rel_nonlin_noncentral' };
21% This defines the number of points used for every algorithm
22indices = { [1:1:17]; [1:1:17]; [1:1:50] };
23% The name of the algorithms in the final plots
24names = { '17pt'; 'nonlin. opt. (17pts)'; 'nonlin. opt. (50pts)' };
25
26% The maximum noise to analyze
27max_noise = 5.0;
28% The step in between different noise levels
29noise_step = 0.1;
30
31%% Run the benchmark
32
33%prepare the overall result arrays
34number_noise_levels = max_noise / noise_step + 1;
35num_algorithms = size(algorithms,1);
36mean_rotation_errors = zeros(num_algorithms,number_noise_levels);
37median_rotation_errors = zeros(num_algorithms,number_noise_levels);
38mean_position_errors = zeros(num_algorithms,number_noise_levels);
39median_position_errors = zeros(num_algorithms,number_noise_levels);
40noise_levels = zeros(1,number_noise_levels);
41
42%Run the experiment
43for n=1:number_noise_levels
44
45    noise = (n - 1) * noise_step;
46    noise_levels(1,n) = noise;
47    display(['Analyzing noise level: ' num2str(noise)])
48
49    rotation_errors = zeros(num_algorithms,iterations);
50    position_errors = zeros(num_algorithms,iterations);
51
52    counter = 0;
53
54    for i=1:iterations
55
56        % generate experiment
57        [v1,v2,t,R] = create2D2DExperiment(pt_number,cam_number,noise,outlier_fraction);
58        [t_perturbed,R_perturbed] = perturb(t,R,0.01);
59        T_perturbed = [R_perturbed,t_perturbed];
60        T_gt = [R,t];
61
62        for a=1:num_algorithms
63            Out = opengv(algorithms{a},indices{a},v1,v2,T_perturbed);
64            [position_error, rotation_error] = evaluateTransformationError( T_gt, Out );
65            position_errors(a,i) = position_error;
66            rotation_errors(a,i) = rotation_error;
67        end
68
69        counter = counter + 1;
70        if counter == 100
71            counter = 0;
72            display(['Iteration ' num2str(i) ' of ' num2str(iterations) '(noise level ' num2str(noise) ')']);
73        end
74    end
75
76    %Now compute the mean and median value of the error for each algorithm
77    for a=1:num_algorithms
78        mean_rotation_errors(a,n) = mean(rotation_errors(a,:));
79        median_rotation_errors(a,n) = median(rotation_errors(a,:));
80        mean_position_errors(a,n) = mean(position_errors(a,:));
81        median_position_errors(a,n) = median(position_errors(a,:));
82    end
83
84end
85
86%% Plot the results
87
88figure(1)
89plot(noise_levels,mean_rotation_errors,'LineWidth',2)
90legend(names,'Location','NorthWest')
91xlabel('noise level [pix]')
92ylabel('mean rot. error [rad]')
93grid on
94
95figure(2)
96plot(noise_levels,median_rotation_errors,'LineWidth',2)
97legend(names,'Location','NorthWest')
98xlabel('noise level [pix]')
99ylabel('median rot. error [rad]')
100grid on
101
102figure(3)
103plot(noise_levels,mean_position_errors,'LineWidth',2)
104legend(names,'Location','NorthWest')
105xlabel('noise level [pix]')
106ylabel('mean pos. error [m]')
107grid on
108
109figure(4)
110plot(noise_levels,median_position_errors,'LineWidth',2)
111legend(names,'Location','NorthWest')
112xlabel('noise level [pix]')
113ylabel('median pos. error [m]')
114grid on