View Javadoc
1   /*
2    * Copyright (C) 2016 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.ar.calibration.estimators;
17  
18  import com.irurueta.ar.calibration.DualAbsoluteQuadric;
19  import com.irurueta.geometry.PinholeCamera;
20  import com.irurueta.geometry.estimators.LockedException;
21  import com.irurueta.geometry.estimators.NotReadyException;
22  import com.irurueta.numerical.robust.RobustEstimatorException;
23  import com.irurueta.numerical.robust.RobustEstimatorMethod;
24  
25  import java.util.List;
26  
27  /**
28   * This is an abstract class for algorithms to robustly find the best
29   * DualAbsoluteQuadric (DAQ) for provided collection of cameras.
30   * Implementations of this class should be able to detect and discard outliers
31   * in order to find the best solution.
32   */
33  @SuppressWarnings("DuplicatedCode")
34  public abstract class DualAbsoluteQuadricRobustEstimator {
35  
36      /**
37       * Default robust estimator method when none is provided.
38       */
39      public static final RobustEstimatorMethod DEFAULT_ROBUST_METHOD = RobustEstimatorMethod.LMEDS;
40  
41      /**
42       * Default amount of progress variation before notifying a change in
43       * estimation progress. By default, this is set to 5%.
44       */
45      public static final float DEFAULT_PROGRESS_DELTA = 0.05f;
46  
47      /**
48       * Minimum allowed value for progress delta.
49       */
50      public static final float MIN_PROGRESS_DELTA = 0.0f;
51  
52      /**
53       * Maximum allowed value for progress delta.
54       */
55      public static final float MAX_PROGRESS_DELTA = 1.0f;
56  
57      /**
58       * Constant defining default confidence of the estimated result, which is
59       * 99%. This means that with a probability of 99% estimation will be
60       * accurate because chosen sub-samples will be inliers.
61       */
62      public static final double DEFAULT_CONFIDENCE = 0.99;
63  
64      /**
65       * Default maximum allowed number of iterations.
66       */
67      public static final int DEFAULT_MAX_ITERATIONS = 5000;
68  
69      /**
70       * Minimum allowed confidence value.
71       */
72      public static final double MIN_CONFIDENCE = 0.0;
73  
74      /**
75       * Maximum allowed confidence value.
76       */
77      public static final double MAX_CONFIDENCE = 1.0;
78  
79      /**
80       * Minimum allowed number of iterations.
81       */
82      public static final int MIN_ITERATIONS = 1;
83  
84      /**
85       * Cameras to estimate dual absolute quadric (DAQ).
86       */
87      protected List<PinholeCamera> cameras;
88  
89      /**
90       * Internal non-robust estimator of DAQ.
91       */
92      protected final LMSEDualAbsoluteQuadricEstimator daqEstimator;
93  
94      /**
95       * Listener to be notified of events such as when estimation starts, ends or
96       * its progress significantly changes.
97       */
98      protected DualAbsoluteQuadricRobustEstimatorListener listener;
99  
100     /**
101      * Indicates if this estimator is locked because an estimation is being
102      * computed.
103      */
104     protected boolean locked;
105 
106     /**
107      * Amount of progress variation before notifying a progress change during
108      * estimation.
109      */
110     protected float progressDelta;
111 
112     /**
113      * Amount of confidence expressed as a value between 0.0 and 1.0 (which is
114      * equivalent to 100%). The amount of confidence indicates the probability
115      * that the estimated result is correct. Usually this value will be close
116      * to 1.0, but not exactly 1.0.
117      */
118     protected double confidence;
119 
120     /**
121      * Maximum allowed number of iterations. When the maximum number of
122      * iterations is exceeded, result will not be available, however an
123      * approximate result will be available for retrieval.
124      */
125     protected int maxIterations;
126 
127     /**
128      * Constructor.
129      */
130     protected DualAbsoluteQuadricRobustEstimator() {
131         progressDelta = DEFAULT_PROGRESS_DELTA;
132         confidence = DEFAULT_CONFIDENCE;
133         maxIterations = DEFAULT_MAX_ITERATIONS;
134         daqEstimator = new LMSEDualAbsoluteQuadricEstimator();
135     }
136 
137     /**
138      * Constructor.
139      *
140      * @param listener listener to be notified of events such as when
141      *                 estimation starts, ends or its progress significantly changes.
142      */
143     protected DualAbsoluteQuadricRobustEstimator(final DualAbsoluteQuadricRobustEstimatorListener listener) {
144         this();
145         this.listener = listener;
146     }
147 
148     /**
149      * Constructor.
150      *
151      * @param cameras list of cameras used to estimate the dual absolute
152      *                quadric (DAQ), which can be used to obtain pinhole camera intrinsic
153      *                parameters.
154      * @throws IllegalArgumentException if not enough cameras are provided
155      *                                  for default settings. Hence, at least 2 cameras must be provided.
156      */
157     protected DualAbsoluteQuadricRobustEstimator(final List<PinholeCamera> cameras) {
158         this();
159         internalSetCameras(cameras);
160     }
161 
162     /**
163      * Constructor.
164      *
165      * @param cameras  list of cameras used to estimate the dual absolute
166      *                 quadric (DAQ), which can be used to obtain pinhole camera intrinsic
167      *                 parameters.
168      * @param listener listener to be notified of events such as when
169      *                 estimation starts, ends or its progress significantly changes.
170      * @throws IllegalArgumentException if not enough cameras are provided
171      *                                  for default settings. Hence, at least 2 cameras must be provided.
172      */
173     protected DualAbsoluteQuadricRobustEstimator(
174             final List<PinholeCamera> cameras, final DualAbsoluteQuadricRobustEstimatorListener listener) {
175         this(listener);
176         internalSetCameras(cameras);
177     }
178 
179     /**
180      * Returns boolean indicating whether camera skewness is assumed to be zero
181      * or not.
182      * Skewness determines whether LCD sensor cells are properly aligned or not,
183      * where zero indicates perfect alignment.
184      * Typically, skewness is a value equal or very close to zero.
185      *
186      * @return true if camera skewness is assumed to be zero, otherwise camera
187      * skewness is estimated.
188      */
189     public boolean isZeroSkewness() {
190         return daqEstimator.isZeroSkewness();
191     }
192 
193     /**
194      * Sets boolean indicating whether camera skewness is assumed to be zero or
195      * not.
196      * Skewness determines whether LCD sensor cells are properly aligned or not,
197      * where zero indicates perfect alignment.
198      * Typically, skewness is a value equal or very close to zero.
199      *
200      * @param zeroSkewness true if camera skewness is assumed to be zero,
201      *                     otherwise camera skewness is estimated.
202      * @throws LockedException if estimator is locked.
203      */
204     public void setZeroSkewness(final boolean zeroSkewness) throws LockedException {
205         if (isLocked()) {
206             throw new LockedException();
207         }
208 
209         daqEstimator.setZeroSkewness(zeroSkewness);
210     }
211 
212     /**
213      * Returns boolean indicating whether principal point is assumed to be at
214      * origin of coordinates or not.
215      * Typically principal point is located at image center (origin of
216      * coordinates), and usually matches the center of radial distortion if
217      * it is taken into account.
218      *
219      * @return true if principal point is assumed to be at origin of
220      * coordinates, false if principal point must be estimated
221      */
222     public boolean isPrincipalPointAtOrigin() {
223         return daqEstimator.isPrincipalPointAtOrigin();
224     }
225 
226     /**
227      * Sets boolean indicating whether principal point is assumed to be at
228      * origin of coordinates or not.
229      * Typically principal point is located at image center (origin of
230      * coordinates), and usually matches the center of radial distortion if it
231      * is taken into account.
232      *
233      * @param principalPointAtOrigin true if principal point is assumed to be at
234      *                               origin of coordinates, false if principal point must be estimated.
235      * @throws LockedException if estimator is locked.
236      */
237     public void setPrincipalPointAtOrigin(final boolean principalPointAtOrigin) throws LockedException {
238         if (isLocked()) {
239             throw new LockedException();
240         }
241 
242         daqEstimator.setPrincipalPointAtOrigin(principalPointAtOrigin);
243     }
244 
245     /**
246      * Returns boolean indicating whether aspect ratio of focal distances (i.e.
247      * vertical focal distance divided by horizontal focal distance) is known or
248      * not.
249      * Notice that focal distance aspect ratio is not related to image size
250      * aspect ratio. Typically, LCD sensor cells are square and hence aspect
251      * ratio of focal distances is known and equal to 1.
252      * This value is only taken into account if skewness is assumed to be zero,
253      * otherwise it is ignored.
254      *
255      * @return true if focal distance aspect ratio is known, false otherwise.
256      */
257     public boolean isFocalDistanceAspectRatioKnown() {
258         return daqEstimator.isFocalDistanceAspectRatioKnown();
259     }
260 
261     /**
262      * Sets value indicating whether aspect ratio of focal distances (i.e.
263      * vertical focal distance divided by horizontal focal distance) is known or
264      * not.
265      * Notice that focal distance aspect ratio is not related to image size
266      * aspect ratio. Typically, LCD sensor cells are square and hence aspect
267      * ratio of focal distances is known and equal to 1.
268      * This value is only taken into account if skewness is assumed to be zero,
269      * otherwise it is ignored.
270      *
271      * @param focalDistanceAspectRatioKnown true if focal distance aspect ratio
272      *                                      is known, false otherwise.
273      * @throws LockedException if estimator is locked.
274      */
275     public void setFocalDistanceAspectRatioKnown(final boolean focalDistanceAspectRatioKnown) throws LockedException {
276         if (isLocked()) {
277             throw new LockedException();
278         }
279 
280         daqEstimator.setFocalDistanceAspectRatioKnown(focalDistanceAspectRatioKnown);
281     }
282 
283     /**
284      * Returns aspect ratio of focal distances (i.e. vertical focal distance
285      * divided by horizontal focal distance).
286      * This value is only taken into account if skewness is assumed to be zero
287      * and focal distance aspect ratio is marked as known, otherwise it is
288      * ignored.
289      * By default, this is 1.0, since it is taken into account that typically
290      * LCD sensor cells are square and hence aspect ratio focal distances is
291      * known and equal to 1.
292      * Notice that focal distance aspect ratio is not related to image size
293      * aspect ratio.
294      * Notice that a negative aspect ratio indicates that vertical axis is
295      * reversed. This can be useful in some situations where image vertical
296      * coordinates are reversed respect to the physical world (i.e. in computer
297      * graphics typically image vertical coordinates go downwards, while in
298      * physical world they go upwards).
299      *
300      * @return aspect ratio of focal distances.
301      */
302     public double getFocalDistanceAspectRatio() {
303         return daqEstimator.getFocalDistanceAspectRatio();
304     }
305 
306     /**
307      * Sets aspect ratio of focal distances (i.e. vertical focal distance
308      * divided by horizontal focal distance).
309      * This value is only taken into account if skewness is assumed to be zero
310      * and focal distance aspect ratio is marked as known, otherwise it is
311      * ignored.
312      * By default, this is 1.0, since it is taken into account that typically
313      * LCD sensor cells are square and hence aspect ratio focal distances is
314      * known and equal to 1.
315      * Notice that focal distance aspect ratio is not related to image size
316      * aspect ratio
317      * Notice that a negative aspect ratio indicates that vertical axis is
318      * reversed. This can be useful in some situations where image vertical
319      * coordinates are reversed respect to the physical world (i.e. in computer
320      * graphics typically image vertical coordinates go downwards, while in
321      * physical world they go upwards).
322      *
323      * @param focalDistanceAspectRatio aspect ratio of focal distances to be set.
324      * @throws LockedException          if estimator is locked.
325      * @throws IllegalArgumentException if focal distance aspect ratio is too
326      *                                  close to zero, as it might produce numerical instabilities.
327      */
328     public void setFocalDistanceAspectRatio(final double focalDistanceAspectRatio) throws LockedException {
329         if (isLocked()) {
330             throw new LockedException();
331         }
332 
333         daqEstimator.setFocalDistanceAspectRatio(focalDistanceAspectRatio);
334     }
335 
336     /**
337      * Indicates whether a singular DAQ is enforced or not.
338      * Dual Absolute Quadric is singular (has rank 3) in any projective space,
339      * however, due to noise in samples, estimated DAQ might not be fully
340      * singular.
341      *
342      * @return true when singular DAQ is enforced, false otherwise.
343      */
344     public boolean isSingularityEnforced() {
345         return daqEstimator.isSingularityEnforced();
346     }
347 
348     /**
349      * Specifies whether a singular DAQ is enforced or not.
350      * Dual Absolute Quadric is singular (has rank 3) in any projective space,
351      * however, due to noise in samples, estimated DAQ might not be fully
352      * singular.
353      *
354      * @param singularityEnforced true when singular DAQ is enforced, false
355      *                            otherwise.
356      * @throws LockedException if estimator is locked.
357      */
358     public void setSingularityEnforced(final boolean singularityEnforced) throws LockedException {
359         if (isLocked()) {
360             throw new LockedException();
361         }
362         daqEstimator.setSingularityEnforced(singularityEnforced);
363     }
364 
365     /**
366      * Indicates whether enforced singularity will be validated by checking that
367      * determinant of estimated Dual Absolute Quadric (DAQ) is below a certain
368      * threshold.
369      *
370      * @return true if enforced singularity is validated, false otherwise.
371      */
372     public boolean isEnforcedSingularityValidated() {
373         return daqEstimator.isEnforcedSingularityValidated();
374     }
375 
376     /**
377      * Specifies whether enforced singularity will be validated by checking that
378      * determinant of estimated Dual Absolute Quadric (DAQ) is below a certain
379      * threshold.
380      *
381      * @param validateEnforcedSingularity true if enforced singularity is
382      *                                    validated, false otherwise.
383      * @throws LockedException if estimator is locked.
384      */
385     public void setEnforcedSingularityValidated(final boolean validateEnforcedSingularity) throws LockedException {
386         if (isLocked()) {
387             throw new LockedException();
388         }
389         daqEstimator.setEnforcedSingularityValidated(validateEnforcedSingularity);
390     }
391 
392     /**
393      * Returns threshold to determine whether estimated Dual Absolute Quadric
394      * (DAQ) has rank 3 or not when validation is enabled.
395      *
396      * @return threshold to determine whether estimated DAQ has rank 3 or not.
397      */
398     public double getDeterminantThreshold() {
399         return daqEstimator.getDeterminantThreshold();
400     }
401 
402     /**
403      * Sets threshold to determine whether estimated Dual Absolute Quadric (DAQ)
404      * has rank 3 or not when validation is enabled.
405      *
406      * @param determinantThreshold threshold to determine whether estimated DAQ
407      *                             has rank 3 or not.
408      * @throws IllegalArgumentException if provided value is zero or negative.
409      * @throws LockedException          if estimator is locked.
410      */
411     public void setDeterminantThreshold(final double determinantThreshold) throws LockedException {
412         if (isLocked()) {
413             throw new LockedException();
414         }
415         daqEstimator.setDeterminantThreshold(determinantThreshold);
416     }
417 
418     /**
419      * Returns reference to listener to be notified of events such as when
420      * estimation starts, ends or its progress significantly changes.
421      *
422      * @return listener to be notified of events.
423      */
424     public DualAbsoluteQuadricRobustEstimatorListener getListener() {
425         return listener;
426     }
427 
428     /**
429      * Sets listener to be notified of events such as when estimation starts,
430      * ends or its progress significantly changes.
431      *
432      * @param listener listener to be notified of events.
433      * @throws LockedException if robust estimator is locked.
434      */
435     public void setListener(final DualAbsoluteQuadricRobustEstimatorListener listener) throws LockedException {
436         if (isLocked()) {
437             throw new LockedException();
438         }
439         this.listener = listener;
440     }
441 
442     /**
443      * Indicates whether listener has been provided and is available for
444      * retrieval.
445      *
446      * @return true if available, false otherwise.
447      */
448     public boolean isListenerAvailable() {
449         return listener != null;
450     }
451 
452     /**
453      * Indicates whether this instance is locked.
454      *
455      * @return true if this estimator is busy estimating the Dual Absolute
456      * Quadric, false otherwise.
457      */
458     public boolean isLocked() {
459         return locked;
460     }
461 
462     /**
463      * Returns amount of progress variation before notifying a progress change
464      * during estimation.
465      *
466      * @return amount of progress variation before notifying a progress change
467      * during estimation.
468      */
469     public float getProgressDelta() {
470         return progressDelta;
471     }
472 
473     /**
474      * Sets amount of progress variation before notifying a progress change
475      * during estimation.
476      *
477      * @param progressDelta amount of progress variation before notifying a
478      *                      progress change during estimation.
479      * @throws IllegalArgumentException if progress delta is less than zero or
480      *                                  greater than 1.
481      * @throws LockedException          if this estimator is locked because an estimation
482      *                                  is being computed.
483      */
484     public void setProgressDelta(final float progressDelta) throws LockedException {
485         if (isLocked()) {
486             throw new LockedException();
487         }
488         if (progressDelta < MIN_PROGRESS_DELTA || progressDelta > MAX_PROGRESS_DELTA) {
489             throw new IllegalArgumentException();
490         }
491         this.progressDelta = progressDelta;
492     }
493 
494     /**
495      * Returns amount of confidence expressed as a value between 0.0 and 1.0
496      * (which is equivalent to 100%). The amount of confidence indicates the
497      * probability that the estimated result is correct. Usually this value will
498      * be close to 1.0, but not exactly 1.0.
499      *
500      * @return amount of confidence as a value between 0.0 and 1.0.
501      */
502     public double getConfidence() {
503         return confidence;
504     }
505 
506     /**
507      * Sets amount of confidence expressed as a value between 0.0 and 1.0 (which
508      * is equivalent to 100%). The amount of confidence indicates the
509      * probability that the estimated result is correct. Usually this value will
510      * be close to 1.0, but not exactly 1.0.
511      *
512      * @param confidence confidence to be set as a value between 0.0 and 1.0.
513      * @throws IllegalArgumentException if provided value is not between 0.0 and
514      *                                  1.0.
515      * @throws LockedException          if this estimator is locked because an estimator
516      *                                  is being computed.
517      */
518     public void setConfidence(final double confidence) throws LockedException {
519         if (isLocked()) {
520             throw new LockedException();
521         }
522         if (confidence < MIN_CONFIDENCE || confidence > MAX_CONFIDENCE) {
523             throw new IllegalArgumentException();
524         }
525         this.confidence = confidence;
526     }
527 
528     /**
529      * Returns maximum allowed number of iterations. If maximum allowed number
530      * of iterations is achieved without converging to a result when calling
531      * estimate(), a RobustEstimatorException will be raised.
532      *
533      * @return maximum allowed number of iterations.
534      */
535     public int getMaxIterations() {
536         return maxIterations;
537     }
538 
539     /**
540      * Sets maximum allowed number of iterations. When the maximum number of
541      * iterations is exceeded, result will not be available, however an
542      * approximate result will be available for retrieval.
543      *
544      * @param maxIterations maximum allowed number of iterations to be set.
545      * @throws IllegalArgumentException if provided value is less than 1.
546      * @throws LockedException          if this estimator is locked because an estimation
547      *                                  is being computed.
548      */
549     public void setMaxIterations(final int maxIterations) throws LockedException {
550         if (isLocked()) {
551             throw new LockedException();
552         }
553         if (maxIterations < MIN_ITERATIONS) {
554             throw new IllegalArgumentException();
555         }
556         this.maxIterations = maxIterations;
557     }
558 
559     /**
560      * Obtains the list of cameras used to estimate the Dual Absolute Quadric
561      * (DAQ).
562      *
563      * @return list of cameras to estimate the DAQ.
564      */
565     public List<PinholeCamera> getCameras() {
566         return cameras;
567     }
568 
569     /**
570      * Sets the list of cameras used to estimate the Dual Absolute Quadric
571      * (DAQ).
572      *
573      * @param cameras list of cameras used to estimate the DAQ.
574      * @throws IllegalArgumentException if list is null.
575      * @throws LockedException          if estimator is locked.
576      */
577     public final void setCameras(final List<PinholeCamera> cameras) throws LockedException {
578         if (isLocked()) {
579             throw new LockedException();
580         }
581         internalSetCameras(cameras);
582     }
583 
584     /**
585      * Returns minimum number of required cameras needed to estimate the
586      * Dual Absolute Quadric (DAQ).
587      * At least 8 equations are needed to solve the DAQ
588      * For each imposed constraint, one less equation is required.
589      * Depending on the number of constraints more or less cameras will be
590      * required.
591      * If zero skewness is enforced, a solution is available with 8 cameras.
592      * If zero skewness and focal distance aspect ratio is known, then a
593      * solution is available with 4 cameras.
594      * If principal point is located at origin, then a solution is available
595      * with 4 cameras.
596      * If zero skewness and principal point at origin are enforced, then a
597      * solution is available with 3 cameras
598      * If zero skewness is enforced, focal distance aspect ratio is known and
599      * principal point is at origin, then a solution is available with 2
600      * cameras.
601      * NOTE: minimum number of cameras considers only the cameras providing
602      * additional information. If a camera is equivalent to another one or does
603      * not provide additional information (such as a camera at the origin with
604      * no rotation), then more cameras will be needed.
605      *
606      * @return minimum number of required cameras needed to estimate the Dual
607      * Absolute Quadric (DAQ) or -1 if constraints configurations is not valid.
608      */
609     public int getMinNumberOfRequiredCameras() {
610         return daqEstimator.getMinNumberOfRequiredCameras();
611     }
612 
613     /**
614      * Indicates whether current constraints are enough to start the estimation.
615      * In order to obtain a linear solution for the DAQ estimation, we need at
616      * least the principal point at origin constraint.
617      *
618      * @return true if constraints are valid, false otherwise.
619      */
620     public boolean areValidConstraints() {
621         return daqEstimator.areValidConstraints();
622     }
623 
624     /**
625      * Returns value indicating whether required data has been provided so that
626      * DAQ estimation can start.
627      * If true, estimator is ready to compute the DAQ, otherwise more data needs
628      * to be provided.
629      *
630      * @return true if estimator is ready, false otherwise.
631      */
632     public boolean isReady() {
633         return cameras != null && cameras.size() >= getMinNumberOfRequiredCameras() && areValidConstraints();
634     }
635 
636     /**
637      * Returns quality scores corresponding to each camera.
638      * The larger the score value the better the quality of the camera.
639      * This implementation always returns null.
640      * Subclasses using quality scores must implement proper behaviour.
641      *
642      * @return quality scores corresponding to each camera.
643      */
644     public double[] getQualityScores() {
645         return null;
646     }
647 
648     /**
649      * Sets quality scores corresponding to each camera.
650      * The larger the score value the better the quality of the camera.
651      * This implementation makes no action.
652      * Subclasses using quality scores must implement proper behaviour.
653      *
654      * @param qualityScores quality scores corresponding to each camera.
655      * @throws LockedException          if robust estimator is locked because an
656      *                                  estimation is already in progress.
657      * @throws IllegalArgumentException if provided quality scores length is
658      *                                  smaller than minimum required number of cameras.
659      */
660     public void setQualityScores(final double[] qualityScores) throws LockedException {
661     }
662 
663     /**
664      * Estimates the Dual Absolute Quadric using provided cameras.
665      *
666      * @return estimated Dual Absolute Quadric (DAQ).
667      * @throws LockedException          if robust estimator is locked.
668      * @throws NotReadyException        if no valid input data has already been
669      *                                  provided.
670      * @throws RobustEstimatorException if estimation fails for any reason
671      *                                  (i.e. numerical instability, no solution available, etc).
672      */
673     public abstract DualAbsoluteQuadric estimate() throws LockedException, NotReadyException, RobustEstimatorException;
674 
675     /**
676      * Returns method being used for robust estimation.
677      *
678      * @return method being used for robust estimation.
679      */
680     public abstract RobustEstimatorMethod getMethod();
681 
682     /**
683      * Creates a dual absolute quadric robust estimator using provided method.
684      *
685      * @param method method of a robust estimator algorithm to estimate best
686      *               DAQ.
687      * @return an instance of a dual absolute quadric robust estimator.
688      */
689     public static DualAbsoluteQuadricRobustEstimator create(final RobustEstimatorMethod method) {
690         return switch (method) {
691             case MSAC -> new MSACDualAbsoluteQuadricRobustEstimator();
692             case RANSAC -> new RANSACDualAbsoluteQuadricRobustEstimator();
693             case PROSAC -> new PROSACDualAbsoluteQuadricRobustEstimator();
694             case PROMEDS -> new PROMedSDualAbsoluteQuadricRobustEstimator();
695             default -> new LMedSDualAbsoluteQuadricRobustEstimator();
696         };
697     }
698 
699     /**
700      * Creates a dual absolute quadric robust estimator using provided
701      * cameras.
702      *
703      * @param cameras       list of cameras.
704      * @param qualityScores quality scores corresponding to each camera.
705      * @param method        method of a robust estimator algorithm to estimate best
706      *                      DAQ.
707      * @return an instance of a dual absolute quadric robust estimator.
708      * @throws IllegalArgumentException if provided list of cameras and quality
709      *                                  scores don't have the same size or size is too short.
710      */
711     public static DualAbsoluteQuadricRobustEstimator create(
712             final List<PinholeCamera> cameras, final double[] qualityScores, final RobustEstimatorMethod method) {
713         return switch (method) {
714             case MSAC -> new MSACDualAbsoluteQuadricRobustEstimator(cameras);
715             case RANSAC -> new RANSACDualAbsoluteQuadricRobustEstimator(cameras);
716             case PROSAC -> new PROSACDualAbsoluteQuadricRobustEstimator(cameras, qualityScores);
717             case PROMEDS -> new PROMedSDualAbsoluteQuadricRobustEstimator(cameras, qualityScores);
718             default -> new LMedSDualAbsoluteQuadricRobustEstimator(cameras);
719         };
720     }
721 
722     /**
723      * Creates a dual absolute quadric robust estimator using provided
724      * cameras.
725      *
726      * @param cameras list of cameras.
727      * @param method  method of a robust estimator algorithm to estimate
728      *                best DAQ.
729      * @return an instance of a dual absolute quadric robust estimator.
730      * @throws IllegalArgumentException if provided list of cameras is too
731      *                                  short.
732      */
733     public static DualAbsoluteQuadricRobustEstimator create(
734             final List<PinholeCamera> cameras, final RobustEstimatorMethod method) {
735         return switch (method) {
736             case MSAC -> new MSACDualAbsoluteQuadricRobustEstimator(cameras);
737             case RANSAC -> new RANSACDualAbsoluteQuadricRobustEstimator(cameras);
738             case PROSAC -> new PROSACDualAbsoluteQuadricRobustEstimator(cameras);
739             case PROMEDS -> new PROMedSDualAbsoluteQuadricRobustEstimator(cameras);
740             default -> new LMedSDualAbsoluteQuadricRobustEstimator(cameras);
741         };
742     }
743 
744     /**
745      * Creates a dual absolute quadric robust estimator using default method.
746      *
747      * @return an instance of a dual absolute quadric robust estimator.
748      */
749     public static DualAbsoluteQuadricRobustEstimator create() {
750         return create(DEFAULT_ROBUST_METHOD);
751     }
752 
753     /**
754      * Creates a dual absolute quadric robust estimator using provided
755      * cameras.
756      *
757      * @param cameras       list of cameras.
758      * @param qualityScores quality scores corresponding to each camera.
759      * @return an instance of a dual absolute quadric robust estimator.
760      * @throws IllegalArgumentException if provided list of cameras and quality
761      *                                  scores don't have the same size or size is too short.
762      */
763     public static DualAbsoluteQuadricRobustEstimator create(
764             final List<PinholeCamera> cameras, final double[] qualityScores) {
765         return create(cameras, qualityScores, DEFAULT_ROBUST_METHOD);
766     }
767 
768     /**
769      * Creates a dual absolute quadric robust estimator using provided
770      * cameras.
771      *
772      * @param cameras list of cameras.
773      * @return an instance of a dual absolute quadric robust estimator.
774      * @throws IllegalArgumentException if provided list of cameras is too
775      *                                  short.
776      */
777     public static DualAbsoluteQuadricRobustEstimator create(final List<PinholeCamera> cameras) {
778         return create(cameras, DEFAULT_ROBUST_METHOD);
779     }
780 
781     /**
782      * Computes the residual between a dual absolute quadric (DAQ) and a pinhole
783      * camera.
784      *
785      * @param daq    a dual absolute quadric (DAQ).
786      * @param camera a camera.
787      * @return residual.
788      */
789     protected double residual(final DualAbsoluteQuadric daq, final PinholeCamera camera) {
790 
791         daq.normalize();
792         camera.normalize();
793         final var cameraMatrix = camera.getInternalMatrix();
794 
795         final var p11 = cameraMatrix.getElementAt(0, 0);
796         final var p21 = cameraMatrix.getElementAt(1, 0);
797         final var p31 = cameraMatrix.getElementAt(2, 0);
798 
799         final var p12 = cameraMatrix.getElementAt(0, 1);
800         final var p22 = cameraMatrix.getElementAt(1, 1);
801         final var p32 = cameraMatrix.getElementAt(2, 1);
802 
803         final var p13 = cameraMatrix.getElementAt(0, 2);
804         final var p23 = cameraMatrix.getElementAt(1, 2);
805         final var p33 = cameraMatrix.getElementAt(2, 2);
806 
807         final var p14 = cameraMatrix.getElementAt(0, 3);
808         final var p24 = cameraMatrix.getElementAt(1, 3);
809         final var p34 = cameraMatrix.getElementAt(2, 3);
810 
811         var residual = 0.0;
812         if (isPrincipalPointAtOrigin()) {
813             if (isZeroSkewness()) {
814                 if (isFocalDistanceAspectRatioKnown()) {
815                     // p2T*daq*p1 = 0
816                     residual += residual2ndRowAnd1stRow(daq, p11, p21, p12, p22, p13, p23, p14, p24);
817 
818                     // p3T*daq*p1 = 0
819                     residual += residual3rdRowAnd1stRow(daq, p11, p31, p12, p32, p13, p33, p14, p34);
820 
821                     // p3T*daw*p2 = 0
822                     residual += residual3rdRowAnd2ndRow(daq, p21, p31, p22, p32, p23, p33, p24, p34);
823 
824                     // p1T*daq*p1*aspectRatio^2 = p2T*daq*p1
825                     residual += residual1stRowEqualTo2ndRow(daq, p11, p21, p12, p22, p13, p23, p14, p24);
826                 } else {
827                     // p2T*daq*p1 = 0
828                     residual += residual2ndRowAnd1stRow(daq, p11, p21, p12, p22, p13, p23, p14, p24);
829 
830                     // p3T*daq*p1 = 0
831                     residual += residual3rdRowAnd1stRow(daq, p11, p31, p12, p32, p13, p33, p14, p34);
832 
833                     //p3T*daq*p2 = 0
834                     residual += residual3rdRowAnd2ndRow(daq, p21, p31, p22, p32, p23, p33, p24, p34);
835                 }
836             } else {
837                 // p3T*daq*p1 = 0
838                 residual += residual3rdRowAnd1stRow(daq, p11, p31, p12, p32, p13, p33, p14, p34);
839 
840                 // p3T*daq*p2 = 0
841                 residual += residual3rdRowAnd2ndRow(daq, p21, p31, p22, p32, p23, p33, p24, p34);
842             }
843         } else {
844             return Double.MAX_VALUE;
845         }
846 
847         return residual;
848     }
849 
850     /**
851      * Computes residual for the equation p2T*daq*p1.
852      *
853      * @param daq estimated DAQ.
854      * @param p11 element (1,1) of camera matrix.
855      * @param p21 element (2,1) of camera matrix.
856      * @param p12 element (1,2) of camera matrix.
857      * @param p22 element (2,2) of camera matrix.
858      * @param p13 element (1,3) of camera matrix.
859      * @param p23 element (2,3) of camera matrix.
860      * @param p14 element (1,4) of camera matrix.
861      * @param p24 element (2,4) of camera matrix.
862      * @return obtained residual (ideally should be zero).
863      */
864     private double residual2ndRowAnd1stRow(
865             final DualAbsoluteQuadric daq, final double p11, final double p21, final double p12, final double p22,
866             final double p13, final double p23, final double p14, final double p24) {
867 
868         final var a = daq.getA();
869         final var b = daq.getB();
870         final var c = daq.getC();
871         final var d = daq.getD();
872         final var e = daq.getE();
873         final var f = daq.getF();
874         final var g = daq.getG();
875         final var h = daq.getH();
876         final var i = daq.getI();
877         final var j = daq.getJ();
878 
879         return a * p21 * p11 + b * p22 * p12 + c * p23 * p13 + d * (p22 * p11 + p21 * p12)
880                 + e * (p23 * p12 + p22 * p13) + f * (p23 * p11 + p21 * p13)
881                 + g * (p24 * p11 + p21 * p14) + h * (p24 * p12 + p22 * p14)
882                 + i * (p24 * p13 + p23 * p14) + j * p24 * p14;
883     }
884 
885     /**
886      * Computes residual for the equation p3T*daq*p1.
887      *
888      * @param daq estimated DAQ.
889      * @param p11 element (1,1) of camera matrix.
890      * @param p31 element (3,1) of camera matrix.
891      * @param p12 element (1,2) of camera matrix.
892      * @param p32 element (3,2) of camera matrix.
893      * @param p13 element (1,3) of camera matrix.
894      * @param p33 element (3,3) of camera matrix.
895      * @param p14 element (1,4) of camera matrix.
896      * @param p34 element (3,4) of camera matrix.
897      * @return obtained residual (ideally should be zero).
898      */
899     private double residual3rdRowAnd1stRow(
900             final DualAbsoluteQuadric daq, final double p11, final double p31, final double p12, final double p32,
901             final double p13, final double p33, final double p14, final double p34) {
902 
903         final var a = daq.getA();
904         final var b = daq.getB();
905         final var c = daq.getC();
906         final var d = daq.getD();
907         final var e = daq.getE();
908         final var f = daq.getF();
909         final var g = daq.getG();
910         final var h = daq.getH();
911         final var i = daq.getI();
912         final var j = daq.getJ();
913 
914         return a * p31 * p11 + b * p32 * p12 + c * p33 * p13 + d * (p32 * p11 + p32 * p12)
915                 + e * (p33 * p12 + p32 * p13) + f * (p33 * p11 + p31 * p13)
916                 + g * (p34 * p11 + p31 * p14) + h * (p34 * p12 + p32 * p14)
917                 + i * (p34 * p13 + p13 * p14) + j * p34 * p14;
918     }
919 
920     /**
921      * Computes residual for the equation p3T*daq*p2.
922      *
923      * @param daq estimated DAQ.
924      * @param p21 element (2,1) of camera matrix.
925      * @param p31 element (3,1) of camera matrix.
926      * @param p22 element (2,2) of camera matrix.
927      * @param p32 element (3,2) of camera matrix.
928      * @param p23 element (2,3) of camera matrix.
929      * @param p33 element (3,3) of camera matrix.
930      * @param p24 element (2,4) of camera matrix.
931      * @param p34 element (3,4) of camera matrix.
932      * @return obtained residual (ideally should be zero).
933      */
934     private double residual3rdRowAnd2ndRow(
935             final DualAbsoluteQuadric daq, final double p21, final double p31, final double p22, final double p32,
936             final double p23, final double p33, final double p24, final double p34) {
937 
938         final var a = daq.getA();
939         final var b = daq.getB();
940         final var c = daq.getC();
941         final var d = daq.getD();
942         final var e = daq.getE();
943         final var f = daq.getF();
944         final var g = daq.getG();
945         final var h = daq.getH();
946         final var i = daq.getI();
947         final var j = daq.getJ();
948 
949         return a * p31 * p21 + b * p32 * p22 + c * p33 * p23 + d * (p32 * p21 + p31 * p22)
950                 + e * (p33 * p22 + p32 * p23) + f * (p33 * p21 + p31 * p23)
951                 + g * (p34 * p21 + p31 * p24) + h * (p34 * p22 + p32 * p24)
952                 + i * (p34 * p23 + p33 * p24) + j * p34 * p24;
953     }
954 
955     /**
956      * Computes residual for the equation p1T*daq*p1 = r^2*p2T*daq*p2.
957      *
958      * @param daq estimated DAQ.
959      * @param p11 element (1,1) of camera matrix.
960      * @param p21 element (2,1) of camera matrix.
961      * @param p12 element (1,2) of camera matrix.
962      * @param p22 element (2,2) of camera matrix.
963      * @param p13 element (1,3) of camera matrix.
964      * @param p23 element (2,3) of camera matrix.
965      * @param p14 element (1,4) of camera matrix.
966      * @param p24 element (2,4) of camera matrix.
967      * @return obtained residual (ideally should be zero).
968      */
969     private double residual1stRowEqualTo2ndRow(
970             final DualAbsoluteQuadric daq, final double p11, final double p21, final double p12, final double p22,
971             final double p13, final double p23, final double p14, final double p24) {
972 
973         final var r = getFocalDistanceAspectRatio();
974         final var r2 = r * r;
975 
976         final var a = daq.getA();
977         final var b = daq.getB();
978         final var c = daq.getC();
979         final var d = daq.getD();
980         final var e = daq.getE();
981         final var f = daq.getF();
982         final var g = daq.getG();
983         final var h = daq.getH();
984         final var i = daq.getI();
985         final var j = daq.getJ();
986 
987         return a * (p11 * p11 * r2 - p21 * p21) + b * (p12 * p12 * r2 - p22 * p22)
988                 + c * (p13 * p13 * r2 - p23 * p23) + d * 2.0 * (p12 * p11 * r2 - p22 * p21)
989                 + e * 2.0 * (p13 * p12 * r2 - p23 * p22) + f * 2.0 * (p13 * p11 * r2 - p23 * p21)
990                 + g * 2.0 * (p14 * p11 * r2 - p24 * p21) + h * 2.0 * (p14 * p12 * r2 - p24 * p22)
991                 + i * 2.0 * (p14 * p13 * r2 - p24 * p23) + j * (p14 * p14 * r2 - p24 * p24);
992     }
993 
994     /**
995      * Sets list of cameras.
996      * This method does not check whether estimator is locked.
997      *
998      * @param cameras list of cameras to estimate DAQ.
999      * @throws IllegalArgumentException if provided list of cameras is null
1000      *                                  or too small.
1001      */
1002     private void internalSetCameras(final List<PinholeCamera> cameras) {
1003         if (cameras == null || cameras.size() < getMinNumberOfRequiredCameras()) {
1004             throw new IllegalArgumentException();
1005         }
1006         this.cameras = cameras;
1007     }
1008 }