View Javadoc
1   /*
2    * Copyright (C) 2020 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.navigation.inertial.calibration.gyroscope;
17  
18  import com.irurueta.navigation.LockedException;
19  import com.irurueta.navigation.NotReadyException;
20  import com.irurueta.navigation.inertial.calibration.CalibrationException;
21  import com.irurueta.navigation.inertial.calibration.StandardDeviationFrameBodyKinematics;
22  import com.irurueta.numerical.robust.LMedSRobustEstimator;
23  import com.irurueta.numerical.robust.LMedSRobustEstimatorListener;
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   * Robustly estimates gyroscope biases, cross couplings and scaling factors
32   * along with G-dependent cross biases introduced on the gyroscope by the
33   * specific forces sensed by the accelerometer using an LMedS algorithm to discard
34   * outliers.
35   * <p>
36   * To use this calibrator at least 7 measurements at different known frames must
37   * be provided. In other words, accelerometer and gyroscope (i.e. body kinematics)
38   * samples must be obtained at 7 different positions, orientations and velocities
39   * (although typically velocities are always zero).
40   * <p>
41   * Measured gyroscope angular rates is assumed to follow the model shown below:
42   * <pre>
43   *     Ωmeas = bg + (I + Mg) * Ωtrue + Gg * ftrue + w
44   * </pre>
45   * Where:
46   * - Ωmeas is the measured gyroscope angular rates. This is a 3x1 vector.
47   * - bg is the gyroscope bias. Ideally, on a perfect gyroscope, this should be a
48   * 3x1 zero vector.
49   * - I is the 3x3 identity matrix.
50   * - Mg is the 3x3 matrix containing cross-couplings and scaling factors. Ideally, on
51   * a perfect gyroscope, this should be a 3x3 zero matrix.
52   * - Ωtrue is ground-truth gyroscope angular rates.
53   * - Gg is the G-dependent cross biases introduced by the specific forces sensed
54   * by the accelerometer. Ideally, on a perfect gyroscope, this should be a 3x3
55   * zero matrix.
56   * - ftrue is ground-truth specific force. This is a 3x1 vector.
57   * - w is measurement noise. This is a 3x1 vector.
58   */
59  public class LMedSRobustKnownFrameGyroscopeCalibrator extends RobustKnownFrameGyroscopeCalibrator {
60  
61      /**
62       * Default value to be used for stop threshold. Stop threshold can be used to
63       * avoid keeping the algorithm unnecessarily iterating in case that best
64       * estimated threshold using median of residuals is not small enough. Once a
65       * solution is found that generates a threshold below this value, the
66       * algorithm will stop.
67       * The stop threshold can be used to prevent the LMedS algorithm iterating
68       * too many times in cases where samples have a very similar accuracy.
69       * For instance, in cases where proportion of outliers is very small (close
70       * to 0%), and samples are very accurate (i.e. 1e-6), the algorithm would
71       * iterate for a long time trying to find the best solution when indeed
72       * there is no need to do that if a reasonable threshold has already been
73       * reached.
74       * Because of this behaviour the stop threshold can be set to a value much
75       * lower than the one typically used in RANSAC, and yet the algorithm could
76       * still produce even smaller thresholds in estimated results.
77       */
78      public static final double DEFAULT_STOP_THRESHOLD = 5e-4;
79  
80      /**
81       * Minimum allowed stop threshold value.
82       */
83      public static final double MIN_STOP_THRESHOLD = 0.0;
84  
85      /**
86       * Threshold to be used to keep the algorithm iterating in case that best
87       * estimated threshold using median of residuals is not small enough. Once
88       * a solution is found that generates a threshold below this value, the
89       * algorithm will stop.
90       * The stop threshold can be used to prevent the LMedS algorithm iterating
91       * too many times in cases where samples have a very similar accuracy.
92       * For instance, in cases where proportion of outliers is very small (close
93       * to 0%), and samples are very accurate (i.e. 1e-6), the algorithm would
94       * iterate for a long time trying to find the best solution when indeed
95       * there is no need to do that if a reasonable threshold has already been
96       * reached.
97       * Because of this behaviour the stop threshold can be set to a value much
98       * lower than the one typically used in RANSAC, and yet the algorithm could
99       * still produce even smaller thresholds in estimated results.
100      */
101     private double stopThreshold = DEFAULT_STOP_THRESHOLD;
102 
103     /**
104      * Constructor.
105      */
106     public LMedSRobustKnownFrameGyroscopeCalibrator() {
107     }
108 
109     /**
110      * Constructor.
111      *
112      * @param listener listener to be notified of events such as when estimation
113      *                 starts, ends or its progress significantly changes.
114      */
115     public LMedSRobustKnownFrameGyroscopeCalibrator(final RobustKnownFrameGyroscopeCalibratorListener listener) {
116         super(listener);
117     }
118 
119     /**
120      * Constructor.
121      *
122      * @param measurements list of body kinematics measurements with standard
123      *                     deviations taken at different frames (positions, orientations
124      *                     and velocities).
125      */
126     public LMedSRobustKnownFrameGyroscopeCalibrator(final List<StandardDeviationFrameBodyKinematics> measurements) {
127         super(measurements);
128     }
129 
130     /**
131      * Constructor.
132      *
133      * @param measurements list of body kinematics measurements with standard
134      *                     deviations taken at different frames (positions, orientations
135      *                     and velocities).
136      * @param listener     listener to be notified of events such as when estimation
137      *                     starts, ends or its progress significantly changes.
138      */
139     public LMedSRobustKnownFrameGyroscopeCalibrator(
140             final List<StandardDeviationFrameBodyKinematics> measurements,
141             final RobustKnownFrameGyroscopeCalibratorListener listener) {
142         super(measurements, listener);
143     }
144 
145     /**
146      * Constructor.
147      *
148      * @param commonAxisUsed indicates whether z-axis is assumed to be common for
149      *                       accelerometer and gyroscope.
150      */
151     public LMedSRobustKnownFrameGyroscopeCalibrator(final boolean commonAxisUsed) {
152         super(commonAxisUsed);
153     }
154 
155     /**
156      * Constructor.
157      *
158      * @param commonAxisUsed indicates whether z-axis is assumed to be common for
159      *                       accelerometer and gyroscope.
160      * @param listener       listener to handle events raised by this calibrator.
161      */
162     public LMedSRobustKnownFrameGyroscopeCalibrator(
163             final boolean commonAxisUsed, final RobustKnownFrameGyroscopeCalibratorListener listener) {
164         super(commonAxisUsed, listener);
165     }
166 
167     /**
168      * Constructor.
169      *
170      * @param measurements   list of body kinematics measurements with standard
171      *                       deviations taken at different frames (positions, orientations
172      *                       and velocities).
173      * @param commonAxisUsed indicates whether z-axis is assumed to be common for
174      *                       accelerometer and gyroscope.
175      */
176     public LMedSRobustKnownFrameGyroscopeCalibrator(
177             final List<StandardDeviationFrameBodyKinematics> measurements, final boolean commonAxisUsed) {
178         super(measurements, commonAxisUsed);
179     }
180 
181     /**
182      * Constructor.
183      *
184      * @param measurements   list of body kinematics measurements with standard
185      *                       deviations taken at different frames (positions, orientations
186      *                       and velocities).
187      * @param commonAxisUsed indicates whether z-axis is assumed to be common for
188      *                       accelerometer and gyroscope.
189      * @param listener       listener to handle events raised by this calibrator.
190      */
191     public LMedSRobustKnownFrameGyroscopeCalibrator(
192             final List<StandardDeviationFrameBodyKinematics> measurements, final boolean commonAxisUsed,
193             final RobustKnownFrameGyroscopeCalibratorListener listener) {
194         super(measurements, commonAxisUsed, listener);
195     }
196 
197     /**
198      * Returns threshold to be used to keep the algorithm iterating in case that
199      * best estimated threshold using median of residuals is not small enough.
200      * Once a solution is found that generates a threshold below this value, the
201      * algorithm will stop.
202      * The stop threshold can be used to prevent the LMedS algorithm to iterate
203      * too many times in cases where samples have a very similar accuracy.
204      * For instance, in cases where proportion of outliers is very small (close
205      * to 0%), and samples are very accurate (i.e. 1e-6), the algorithm would
206      * iterate for a long time trying to find the best solution when indeed
207      * there is no need to do that if a reasonable threshold has already been
208      * reached.
209      * Because of this behaviour the stop threshold can be set to a value much
210      * lower than the one typically used in RANSAC, and yet the algorithm could
211      * still produce even smaller thresholds in estimated results.
212      *
213      * @return stop threshold to stop the algorithm prematurely when a certain
214      * accuracy has been reached.
215      */
216     public double getStopThreshold() {
217         return stopThreshold;
218     }
219 
220     /**
221      * Sets threshold to be used to keep the algorithm iterating in case that
222      * best estimated threshold using median of residuals is not small enough.
223      * Once a solution is found that generates a threshold below this value,
224      * the algorithm will stop.
225      * The stop threshold can be used to prevent the LMedS algorithm to iterate
226      * too many times in cases where samples have a very similar accuracy.
227      * For instance, in cases where proportion of outliers is very small (close
228      * to 0%), and samples are very accurate (i.e. 1e-6), the algorithm would
229      * iterate for a long time trying to find the best solution when indeed
230      * there is no need to do that if a reasonable threshold has already been
231      * reached.
232      * Because of this behaviour the stop threshold can be set to a value much
233      * lower than the one typically used in RANSAC, and yet the algorithm could
234      * still produce even smaller thresholds in estimated results.
235      *
236      * @param stopThreshold stop threshold to stop the algorithm prematurely
237      *                      when a certain accuracy has been reached.
238      * @throws IllegalArgumentException if provided value is zero or negative.
239      * @throws LockedException          if calibrator is currently running.
240      */
241     public void setStopThreshold(final double stopThreshold) throws LockedException {
242         if (running) {
243             throw new LockedException();
244         }
245         if (stopThreshold <= MIN_STOP_THRESHOLD) {
246             throw new IllegalArgumentException();
247         }
248 
249         this.stopThreshold = stopThreshold;
250     }
251 
252     /**
253      * Estimates accelerometer calibration parameters containing bias, scale factors
254      * and cross-coupling errors.
255      *
256      * @throws LockedException      if calibrator is currently running.
257      * @throws NotReadyException    if calibrator is not ready.
258      * @throws CalibrationException if estimation fails for numerical reasons.
259      */
260     @SuppressWarnings("DuplicatedCode")
261     @Override
262     public void calibrate() throws LockedException, NotReadyException, CalibrationException {
263         if (running) {
264             throw new LockedException();
265         }
266         if (!isReady()) {
267             throw new NotReadyException();
268         }
269 
270         final var innerEstimator = new LMedSRobustEstimator<>(new LMedSRobustEstimatorListener<PreliminaryResult>() {
271             @Override
272             public int getTotalSamples() {
273                 return measurements.size();
274             }
275 
276             @Override
277             public int getSubsetSize() {
278                 return preliminarySubsetSize;
279             }
280 
281             @Override
282             public void estimatePreliminarSolutions(
283                     final int[] samplesIndices, final List<PreliminaryResult> solutions) {
284                 computePreliminarySolutions(samplesIndices, solutions);
285             }
286 
287             @Override
288             public double computeResidual(final PreliminaryResult currentEstimation, final int i) {
289                 return computeError(measurements.get(i), currentEstimation);
290             }
291 
292             @Override
293             public boolean isReady() {
294                 return LMedSRobustKnownFrameGyroscopeCalibrator.super.isReady();
295             }
296 
297             @Override
298             public void onEstimateStart(final RobustEstimator<PreliminaryResult> estimator) {
299                 // no action needed
300             }
301 
302             @Override
303             public void onEstimateEnd(final RobustEstimator<PreliminaryResult> estimator) {
304                 // no action needed
305             }
306 
307             @Override
308             public void onEstimateNextIteration(
309                     final RobustEstimator<PreliminaryResult> estimator, final int iteration) {
310                 if (listener != null) {
311                     listener.onCalibrateNextIteration(
312                             LMedSRobustKnownFrameGyroscopeCalibrator.this, iteration);
313                 }
314             }
315 
316             @Override
317             public void onEstimateProgressChange(
318                     final RobustEstimator<PreliminaryResult> estimator, final float progress) {
319                 if (listener != null) {
320                     listener.onCalibrateProgressChange(
321                             LMedSRobustKnownFrameGyroscopeCalibrator.this, progress);
322                 }
323             }
324         });
325 
326         try {
327             running = true;
328 
329             if (listener != null) {
330                 listener.onCalibrateStart(this);
331             }
332 
333             inliersData = null;
334             innerEstimator.setConfidence(confidence);
335             innerEstimator.setMaxIterations(maxIterations);
336             innerEstimator.setProgressDelta(progressDelta);
337             innerEstimator.setStopThreshold(stopThreshold);
338             final var preliminaryResult = innerEstimator.estimate();
339             inliersData = innerEstimator.getInliersData();
340 
341             attemptRefine(preliminaryResult);
342 
343             if (listener != null) {
344                 listener.onCalibrateEnd(this);
345             }
346 
347         } catch (final com.irurueta.numerical.LockedException e) {
348             throw new LockedException(e);
349         } catch (final com.irurueta.numerical.NotReadyException e) {
350             throw new NotReadyException(e);
351         } catch (final RobustEstimatorException e) {
352             throw new CalibrationException(e);
353         } finally {
354             running = false;
355         }
356     }
357 
358     /**
359      * Returns method being used for robust estimation.
360      *
361      * @return method being used for robust estimation.
362      */
363     @Override
364     public RobustEstimatorMethod getMethod() {
365         return RobustEstimatorMethod.LMEDS;
366     }
367 
368     /**
369      * Indicates whether this calibrator requires quality scores for each
370      * measurement/sequence or not.
371      *
372      * @return true if quality scores are required, false otherwise.
373      */
374     @Override
375     public boolean isQualityScoresRequired() {
376         return false;
377     }
378 }