1 /*
2 * Copyright (C) 2015 Alberto Irurueta Carro (alberto@irurueta.com)
3 *
4 * Licensed under the Apache License, Version 2.0 (the "License");
5 * you may not use this file except in compliance with the License.
6 * You may obtain a copy of the License at
7 *
8 * http://www.apache.org/licenses/LICENSE-2.0
9 *
10 * Unless required by applicable law or agreed to in writing, software
11 * distributed under the License is distributed on an "AS IS" BASIS,
12 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13 * See the License for the specific language governing permissions and
14 * limitations under the License.
15 */
16 package com.irurueta.geometry.estimators;
17
18 import com.irurueta.geometry.AffineTransformation3D;
19 import com.irurueta.geometry.CoincidentPointsException;
20 import com.irurueta.geometry.CoordinatesType;
21 import com.irurueta.geometry.Point3D;
22 import com.irurueta.numerical.robust.PROMedSRobustEstimator;
23 import com.irurueta.numerical.robust.PROMedSRobustEstimatorListener;
24 import com.irurueta.numerical.robust.RobustEstimator;
25 import com.irurueta.numerical.robust.RobustEstimatorException;
26 import com.irurueta.numerical.robust.RobustEstimatorMethod;
27
28 import java.util.List;
29
30 /**
31 * Finds the best affine 3D transformation for provided collections of matched
32 * 3D points using PROSAC algorithm.
33 */
34 public class PROMedSPointCorrespondenceAffineTransformation3DRobustEstimator
35 extends PointCorrespondenceAffineTransformation3DRobustEstimator {
36
37
38 /**
39 * Default value to be used for stop threshold. Stop threshold can be used
40 * to keep the algorithm iterating in case that best estimated threshold
41 * using median of residuals is not small enough. Once a solution is found
42 * that generates a threshold below this value, the algorithm will stop.
43 * The stop threshold can be used to prevent the LMedS algorithm iterating
44 * too many times in cases where samples have a very similar accuracy.
45 * For instance, in cases where proportion of outliers is very small (close
46 * to 0%), and samples are very accurate (i.e. 1e-6), the algorithm would
47 * iterate for a long time trying to find the best solution when indeed
48 * there is no need to do that if a reasonable threshold has already been
49 * reached.
50 * Because of this behaviour the stop threshold can be set to a value much
51 * lower than the one typically used in RANSAC, and yet the algorithm could
52 * still produce even smaller thresholds in estimated results.
53 */
54 public static final double DEFAULT_STOP_THRESHOLD = 1.0;
55
56 /**
57 * Minimum allowed stop threshold value.
58 */
59 public static final double MIN_STOP_THRESHOLD = 0.0;
60
61 /**
62 * Threshold to be used to keep the algorithm iterating in case that best
63 * estimated threshold using median of residuals is not small enough. Once
64 * a solution is found that generates a threshold below this value, the
65 * algorithm will stop.
66 * The stop threshold can be used to prevent the LMedS algorithm iterating
67 * too many times in cases where samples have a very similar accuracy.
68 * For instance, in cases where proportion of outliers is very small (close
69 * to 0%), and samples are very accurate (i.e. 1e-6), the algorithm would
70 * iterate for a long time trying to find the best solution when indeed
71 * there is no need to do that if a reasonable threshold has already been
72 * reached.
73 * Because of this behaviour the stop threshold can be set to a value much
74 * lower than the one typically used in RANSAC, and yet the algorithm could
75 * still produce even smaller thresholds in estimated results.
76 */
77 private double stopThreshold;
78
79 /**
80 * Quality scores corresponding to each pair of matched points.
81 * The larger the score value the better the quality of the matching.
82 */
83 private double[] qualityScores;
84
85 /**
86 * Constructor.
87 */
88 public PROMedSPointCorrespondenceAffineTransformation3DRobustEstimator() {
89 super();
90 stopThreshold = DEFAULT_STOP_THRESHOLD;
91 }
92
93 /**
94 * Constructor with lists of points to be used to estimate an affine 3D
95 * transformation.
96 * Points in the list located at the same position are considered to be
97 * matched. Hence, both lists must have the same size, and their size must
98 * be greater or equal than MINIMUM_SIZE.
99 *
100 * @param inputPoints list of input points to be used to estimate an
101 * affine 3D transformation.
102 * @param outputPoints list of output points to be used to estimate an
103 * affine 3D transformation.
104 * @throws IllegalArgumentException if provided lists of points don't have
105 * the same size or their size is smaller than MINIMUM_SIZE.
106 */
107 public PROMedSPointCorrespondenceAffineTransformation3DRobustEstimator(
108 final List<Point3D> inputPoints, final List<Point3D> outputPoints) {
109 super(inputPoints, outputPoints);
110 stopThreshold = DEFAULT_STOP_THRESHOLD;
111 }
112
113 /**
114 * Constructor.
115 *
116 * @param listener listener to be notified of events such as when estimation
117 * starts, ends or its progress significantly changes.
118 */
119 public PROMedSPointCorrespondenceAffineTransformation3DRobustEstimator(
120 final AffineTransformation3DRobustEstimatorListener listener) {
121 super(listener);
122 stopThreshold = DEFAULT_STOP_THRESHOLD;
123 }
124
125 /**
126 * Constructor with listener and lists of points to be used to estimate an
127 * affine 3D transformation.
128 * Points in the list located at the same position are considered to be
129 * matched. Hence, both lists must have the same size, and their size must
130 * be greater or equal than MINIMUM_SIZE.
131 *
132 * @param listener listener to be notified of events such as when estimation
133 * stars, ends or its progress significantly changes.
134 * @param inputPoints list of input points to be used to estimate an
135 * affine 3D transformation.
136 * @param outputPoints list of output points to be used to estimate an
137 * affine 3D transformation.
138 * @throws IllegalArgumentException if provided lists of points don't have
139 * the same size or their size is smaller than MINIMUM_SIZE.
140 */
141 public PROMedSPointCorrespondenceAffineTransformation3DRobustEstimator(
142 final AffineTransformation3DRobustEstimatorListener listener,
143 final List<Point3D> inputPoints, final List<Point3D> outputPoints) {
144 super(listener, inputPoints, outputPoints);
145 stopThreshold = DEFAULT_STOP_THRESHOLD;
146 }
147
148 /**
149 * Constructor.
150 *
151 * @param qualityScores quality scores corresponding to each pair of matched
152 * points.
153 * @throws IllegalArgumentException if provided quality scores length is
154 * smaller than MINIMUM_SIZE (i.e. 3 samples).
155 */
156 public PROMedSPointCorrespondenceAffineTransformation3DRobustEstimator(final double[] qualityScores) {
157 super();
158 stopThreshold = DEFAULT_STOP_THRESHOLD;
159 internalSetQualityScores(qualityScores);
160 }
161
162 /**
163 * Constructor with lists of points to be used to estimate an affine 3D
164 * transformation.
165 * Points in the list located at the same position are considered to be
166 * matched. Hence, both lists must have the same size, and their size must
167 * be greater or equal than MINIMUM_SIZE.
168 *
169 * @param inputPoints list of input points to be used to estimate an
170 * affine 3D transformation.
171 * @param outputPoints list of output points to be used to estimate an
172 * affine 3D transformation.
173 * @param qualityScores quality scores corresponding to each pair of matched
174 * points.
175 * @throws IllegalArgumentException if provided lists of points and array
176 * of quality scores don't have the same size or their size is smaller than
177 * MINIMUM_SIZE.
178 */
179 public PROMedSPointCorrespondenceAffineTransformation3DRobustEstimator(
180 final List<Point3D> inputPoints, final List<Point3D> outputPoints, final double[] qualityScores) {
181 super(inputPoints, outputPoints);
182
183 if (qualityScores.length != inputPoints.size()) {
184 throw new IllegalArgumentException();
185 }
186
187 stopThreshold = DEFAULT_STOP_THRESHOLD;
188 internalSetQualityScores(qualityScores);
189 }
190
191 /**
192 * Constructor.
193 *
194 * @param listener listener to be notified of events such as when estimation
195 * starts, ends or its progress significantly changes
196 * @param qualityScores quality scores corresponding to each pair of matched
197 * points.
198 * @throws IllegalArgumentException if provided quality scores length is
199 * smaller than MINIMUM_SIZE (i.e. 3 samples).
200 */
201 public PROMedSPointCorrespondenceAffineTransformation3DRobustEstimator(
202 final AffineTransformation3DRobustEstimatorListener listener, final double[] qualityScores) {
203 super(listener);
204 stopThreshold = DEFAULT_STOP_THRESHOLD;
205 internalSetQualityScores(qualityScores);
206 }
207
208 /**
209 * Constructor with listener and lists of points to be used to estimate an
210 * affine 3D transformation.
211 * Points in the list located at the same position are considered to be
212 * matched. Hence, both lists must have the same size, and their size must
213 * be greater or equal than MINIMUM_SIZE.
214 *
215 * @param listener listener to be notified of events such as when estimation
216 * stars, ends or its progress significantly changes.
217 * @param inputPoints list of input points to be used to estimate an
218 * affine 3D transformation.
219 * @param outputPoints list of output points to be used to estimate an
220 * affine 3D transformation.
221 * @param qualityScores quality scores corresponding to each pair of matched
222 * points.
223 * @throws IllegalArgumentException if provided lists of points don't have
224 * the same size or their size is smaller than MINIMUM_SIZE.
225 */
226 public PROMedSPointCorrespondenceAffineTransformation3DRobustEstimator(
227 final AffineTransformation3DRobustEstimatorListener listener, final List<Point3D> inputPoints,
228 final List<Point3D> outputPoints, final double[] qualityScores) {
229 super(listener, inputPoints, outputPoints);
230
231 if (qualityScores.length != inputPoints.size()) {
232 throw new IllegalArgumentException();
233 }
234
235 stopThreshold = DEFAULT_STOP_THRESHOLD;
236 internalSetQualityScores(qualityScores);
237 }
238
239 /**
240 * Returns threshold to be used to keep the algorithm iterating in case that
241 * best estimated threshold using median of residuals is not small enough.
242 * Once a solution is found that generates a threshold below this value, the
243 * algorithm will stop.
244 * As in LMedS, the stop threshold can be used to prevent the PROMedS
245 * algorithm iterating too many times in cases where samples have a very
246 * similar accuracy.
247 * For instance, in cases where proportion of outliers is very small (close
248 * to 0%), and samples are very accurate (i.e. 1e-6), the algorithm would
249 * iterate for a long time trying to find the best solution when indeed
250 * there is no need to do that if a reasonable threshold has already been
251 * reached.
252 * Because of this behaviour the stop threshold can be set to a value much
253 * lower than the one typically used in RANSAC, and yet the algorithm could
254 * still produce even smaller thresholds in estimated results.
255 *
256 * @return stop threshold to stop the algorithm prematurely when a certain
257 * accuracy has been reached.
258 */
259 public double getStopThreshold() {
260 return stopThreshold;
261 }
262
263 /**
264 * Sets threshold to be used to keep the algorithm iterating in case that
265 * best estimated threshold using median of residuals is not small enough.
266 * Once a solution is found that generates a threshold below this value, the
267 * algorithm will stop.
268 * As in LMedS, the stop threshold can be used to prevent the PROMedS
269 * algorithm iterating too many times in cases where samples have a very
270 * similar accuracy.
271 * For instance, in cases where proportion of outliers is very small (close
272 * to 0%), and samples are very accurate (i.e. 1e-6), the algorithm would
273 * iterate for a long time trying to find the best solution when indeed
274 * there is no need to do that if a reasonable threshold has already been
275 * reached.
276 * Because of this behaviour the stop threshold can be set to a value much
277 * lower than the one typically used in RANSAC, and yet the algorithm could
278 * still produce even smaller thresholds in estimated results.
279 *
280 * @param stopThreshold stop threshold to stop the algorithm prematurely
281 * when a certain accuracy has been reached.
282 * @throws IllegalArgumentException if provided value is zero or negative.
283 * @throws LockedException if robust estimator is locked because an
284 * estimation is already in progress.
285 */
286 public void setStopThreshold(final double stopThreshold) throws LockedException {
287 if (isLocked()) {
288 throw new LockedException();
289 }
290 if (stopThreshold <= MIN_STOP_THRESHOLD) {
291 throw new IllegalArgumentException();
292 }
293
294 this.stopThreshold = stopThreshold;
295 }
296
297 /**
298 * Returns quality scores corresponding to each pair of matched points.
299 * The larger the score value the better the quality of the matching.
300 *
301 * @return quality scores corresponding to each pair of matched points.
302 */
303 @Override
304 public double[] getQualityScores() {
305 return qualityScores;
306 }
307
308 /**
309 * Sets quality scores corresponding to each pair of matched points.
310 * The larger the score value the better the quality of the matching.
311 *
312 * @param qualityScores quality scores corresponding to each pair of matched
313 * points.
314 * @throws LockedException if robust estimator is locked because an
315 * estimation is already in progress.
316 * @throws IllegalArgumentException if provided quality scores length is
317 * smaller than MINIMUM_SIZE (i.e. 3 samples).
318 */
319 @Override
320 public void setQualityScores(final double[] qualityScores) throws LockedException {
321 if (isLocked()) {
322 throw new LockedException();
323 }
324 internalSetQualityScores(qualityScores);
325 }
326
327 /**
328 * Indicates if estimator is ready to start the affine 2D transformation
329 * estimation.
330 * This is true when input data (i.e. lists of matched points and quality
331 * scores) are provided and a minimum of MINIMUM_SIZE points are available.
332 *
333 * @return true if estimator is ready, false otherwise.
334 */
335 @Override
336 public boolean isReady() {
337 return super.isReady() && qualityScores != null && qualityScores.length == inputPoints.size();
338 }
339
340 /**
341 * Estimates an affine 2D transformation using a robust estimator and
342 * the best set of matched 2D point correspondences found using the robust
343 * estimator.
344 *
345 * @return an affine 2D transformation.
346 * @throws LockedException if robust estimator is locked because an
347 * estimation is already in progress.
348 * @throws NotReadyException if provided input data is not enough to start
349 * the estimation.
350 * @throws RobustEstimatorException if estimation fails for any reason
351 * (i.e. numerical instability, no solution available, etc).
352 */
353 @SuppressWarnings("DuplicatedCode")
354 @Override
355 public AffineTransformation3D estimate() throws LockedException, NotReadyException, RobustEstimatorException {
356 if (isLocked()) {
357 throw new LockedException();
358 }
359 if (!isReady()) {
360 throw new NotReadyException();
361 }
362
363 final var innerEstimator = new PROMedSRobustEstimator<>(
364 new PROMedSRobustEstimatorListener<AffineTransformation3D>() {
365
366 // point to be reused when computing residuals
367 private final Point3D testPoint = Point3D.create(CoordinatesType.HOMOGENEOUS_COORDINATES);
368
369 @Override
370 public double getThreshold() {
371 return stopThreshold;
372 }
373
374 @Override
375 public int getTotalSamples() {
376 return inputPoints.size();
377 }
378
379 @Override
380 public int getSubsetSize() {
381 return AffineTransformation3DRobustEstimator.MINIMUM_SIZE;
382 }
383
384 @Override
385 public void estimatePreliminarSolutions(
386 final int[] samplesIndices, final List<AffineTransformation3D> solutions) {
387 final var inputPoint1 = inputPoints.get(samplesIndices[0]);
388 final var inputPoint2 = inputPoints.get(samplesIndices[1]);
389 final var inputPoint3 = inputPoints.get(samplesIndices[2]);
390 final var inputPoint4 = inputPoints.get(samplesIndices[3]);
391
392 final var outputPoint1 = outputPoints.get(samplesIndices[0]);
393 final var outputPoint2 = outputPoints.get(samplesIndices[1]);
394 final var outputPoint3 = outputPoints.get(samplesIndices[2]);
395 final var outputPoint4 = outputPoints.get(samplesIndices[3]);
396
397 try {
398 final var transformation = new AffineTransformation3D(inputPoint1, inputPoint2, inputPoint3,
399 inputPoint4, outputPoint1, outputPoint2, outputPoint3, outputPoint4);
400 solutions.add(transformation);
401 } catch (final CoincidentPointsException e) {
402 // if points are coincident, no solution is added
403 }
404 }
405
406 @Override
407 public double computeResidual(final AffineTransformation3D currentEstimation, final int i) {
408 final var inputPoint = inputPoints.get(i);
409 final var outputPoint = outputPoints.get(i);
410
411 // transform input point and store result in mTestPoint
412 currentEstimation.transform(inputPoint, testPoint);
413
414 return outputPoint.distanceTo(testPoint);
415 }
416
417 @Override
418 public boolean isReady() {
419 return PROMedSPointCorrespondenceAffineTransformation3DRobustEstimator.this.isReady();
420 }
421
422 @Override
423 public void onEstimateStart(final RobustEstimator<AffineTransformation3D> estimator) {
424 if (listener != null) {
425 listener.onEstimateStart(
426 PROMedSPointCorrespondenceAffineTransformation3DRobustEstimator.this);
427 }
428 }
429
430 @Override
431 public void onEstimateEnd(final RobustEstimator<AffineTransformation3D> estimator) {
432 if (listener != null) {
433 listener.onEstimateEnd(
434 PROMedSPointCorrespondenceAffineTransformation3DRobustEstimator.this);
435 }
436 }
437
438 @Override
439 public void onEstimateNextIteration(
440 final RobustEstimator<AffineTransformation3D> estimator, final int iteration) {
441 if (listener != null) {
442 listener.onEstimateNextIteration(
443 PROMedSPointCorrespondenceAffineTransformation3DRobustEstimator.this,
444 iteration);
445 }
446 }
447
448 @Override
449 public void onEstimateProgressChange(
450 final RobustEstimator<AffineTransformation3D> estimator, final float progress) {
451 if (listener != null) {
452 listener.onEstimateProgressChange(
453 PROMedSPointCorrespondenceAffineTransformation3DRobustEstimator.this,
454 progress);
455 }
456 }
457
458 @Override
459 public double[] getQualityScores() {
460 return qualityScores;
461 }
462 });
463
464 try {
465 locked = true;
466 inliersData = null;
467 innerEstimator.setConfidence(confidence);
468 innerEstimator.setMaxIterations(maxIterations);
469 innerEstimator.setProgressDelta(progressDelta);
470 final var transformation = innerEstimator.estimate();
471 inliersData = innerEstimator.getInliersData();
472 return attemptRefine(transformation);
473 } catch (final com.irurueta.numerical.LockedException e) {
474 throw new LockedException(e);
475 } catch (final com.irurueta.numerical.NotReadyException e) {
476 throw new NotReadyException(e);
477 } finally {
478 locked = false;
479 }
480 }
481
482 /**
483 * Returns method being used for robust estimation.
484 *
485 * @return method being used for robust estimation.
486 */
487 @Override
488 public RobustEstimatorMethod getMethod() {
489 return RobustEstimatorMethod.PROMEDS;
490 }
491
492 /**
493 * Gets standard deviation used for Levenberg-Marquardt fitting during
494 * refinement.
495 * Returned value gives an indication of how much variance each residual
496 * has.
497 * Typically, this value is related to the threshold used on each robust
498 * estimation, since residuals of found inliers are within the range of such
499 * threshold.
500 *
501 * @return standard deviation used for refinement.
502 */
503 @Override
504 protected double getRefinementStandardDeviation() {
505 final var inliersData = (PROMedSRobustEstimator.PROMedSInliersData) getInliersData();
506
507 // avoid setting a threshold too strict
508 final var threshold = inliersData.getEstimatedThreshold();
509 return Math.max(threshold, stopThreshold);
510 }
511
512 /**
513 * Sets quality scores corresponding to each pair of matched points.
514 * This method is used internally and does not check whether instance is
515 * locked or not.
516 *
517 * @param qualityScores quality scores to be set.
518 * @throws IllegalArgumentException if provided quality scores length is
519 * smaller than MINIMUM_SIZE.
520 */
521 private void internalSetQualityScores(final double[] qualityScores) {
522 if (qualityScores.length < MINIMUM_SIZE) {
523 throw new IllegalArgumentException();
524 }
525
526 this.qualityScores = qualityScores;
527 }
528 }