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.algebra.AlgebraException;
19 import com.irurueta.algebra.ArrayUtils;
20 import com.irurueta.algebra.Matrix;
21 import com.irurueta.algebra.Utils;
22 import com.irurueta.algebra.WrongSizeException;
23 import com.irurueta.geometry.InhomogeneousPoint3D;
24 import com.irurueta.geometry.Quaternion;
25 import com.irurueta.geometry.RotationException;
26 import com.irurueta.navigation.LockedException;
27 import com.irurueta.navigation.NotReadyException;
28 import com.irurueta.navigation.inertial.BodyKinematics;
29 import com.irurueta.navigation.inertial.INSLooselyCoupledKalmanInitializerConfig;
30 import com.irurueta.navigation.inertial.INSTightlyCoupledKalmanInitializerConfig;
31 import com.irurueta.navigation.inertial.calibration.AccelerationFixer;
32 import com.irurueta.navigation.inertial.calibration.AngularSpeedTriad;
33 import com.irurueta.navigation.inertial.calibration.BodyKinematicsSequence;
34 import com.irurueta.navigation.inertial.calibration.CalibrationException;
35 import com.irurueta.navigation.inertial.calibration.GyroscopeBiasUncertaintySource;
36 import com.irurueta.navigation.inertial.calibration.GyroscopeCalibrationSource;
37 import com.irurueta.navigation.inertial.calibration.StandardDeviationTimedBodyKinematics;
38 import com.irurueta.numerical.EvaluationException;
39 import com.irurueta.numerical.GradientEstimator;
40 import com.irurueta.numerical.fitting.FittingException;
41 import com.irurueta.numerical.fitting.LevenbergMarquardtMultiDimensionFitter;
42 import com.irurueta.numerical.fitting.LevenbergMarquardtMultiDimensionFunctionEvaluator;
43 import com.irurueta.statistics.MaxIterationsExceededException;
44 import com.irurueta.units.Acceleration;
45 import com.irurueta.units.AccelerationConverter;
46 import com.irurueta.units.AccelerationUnit;
47 import com.irurueta.units.AngularSpeed;
48 import com.irurueta.units.AngularSpeedConverter;
49 import com.irurueta.units.AngularSpeedUnit;
50
51 import java.util.ArrayList;
52 import java.util.List;
53
54 /**
55 * Estimates gyroscope biases, cross couplings, scaling factors and
56 * G-dependent cross biases introduced on the gyroscope by the specific
57 * forces sensed by the accelerometer.
58 * <p>
59 * This calibrator assumes that the IMU is at a more or less fixed location on
60 * Earth, and evaluates sequences of measured body kinematics to perform
61 * calibration for unknown orientations on those provided sequences.
62 * Each provided sequence will be preceded by a static period where mean
63 * specific force will be measured to determine gravity (and hence partial
64 * body attitude).
65 * <p>
66 * To use this calibrator at least 10 sequences are needed (each containing
67 * at least 3 body kinematics measurements) when common z-axis is assumed and
68 * G-dependant cross biases are ignored, otherwise at least 13 sequences are
69 * required (each containing at least 3 body kinematics measurements) when
70 * common z-axis is not assumed.
71 * If G-dependent cross biases are being estimated, then at least 19
72 * measurements are needed when common z-axis is assumed, otherwise at least
73 * 22 sequences are required (each containing at least 3 body kinematics
74 * measurements) when common z-axis is not assumed.
75 * <p>
76 * Measured gyroscope angular rates is assumed to follow the model shown below:
77 * <pre>
78 * Ωmeas = bg + (I + Mg) * Ωtrue + Gg * ftrue + w
79 * </pre>
80 * Where:
81 * - Ωmeas is the measured gyroscope angular rates. This is a 3x1 vector.
82 * - bg is the gyroscope bias. Ideally, on a perfect gyroscope, this should be a
83 * 3x1 zero vector.
84 * - I is the 3x3 identity matrix.
85 * - Mg is the 3x3 matrix containing cross-couplings and scaling factors. Ideally, on
86 * a perfect gyroscope, this should be a 3x3 zero matrix.
87 * - Ωtrue is ground-truth gyroscope angular rates.
88 * - Gg is the G-dependent cross biases introduced by the specific forces sensed
89 * by the accelerometer. Ideally, on a perfect gyroscope, this should be a 3x3
90 * zero matrix.
91 * - ftrue is ground-truth specific force. This is a 3x1 vector.
92 * - w is measurement noise. This is a 3x1 vector.
93 * <p>
94 * This calibrator is based on the ideas of
95 * David Tedaldi, Alberto Pretto, Emmanuelle Menegatti. A Robust and Easy to
96 * Implement Method for IMU Calibration without External Equipments.
97 */
98 public class EasyGyroscopeCalibrator implements
99 GyroscopeNonLinearCalibrator, UnknownBiasGyroscopeCalibrator, GyroscopeCalibrationSource,
100 GyroscopeBiasUncertaintySource, OrderedBodyKinematicsSequenceGyroscopeCalibrator,
101 AccelerometerDependentGyroscopeCalibrator {
102
103 /**
104 * Indicates whether by default a common z-axis is assumed for both the accelerometer
105 * and gyroscope.
106 */
107 public static final boolean DEFAULT_USE_COMMON_Z_AXIS = true;
108
109 /**
110 * Indicates that by default G-dependent cross biases introduced
111 * by the accelerometer on the gyroscope are estimated.
112 */
113 public static final boolean DEFAULT_ESTIMATE_G_DEPENDENT_CROSS_BIASES = true;
114
115 /**
116 * Number of unknowns when common z-axis is assumed for both the accelerometer
117 * and gyroscope when G-dependent cross biases are being estimated.
118 */
119 public static final int COMMON_Z_AXIS_UNKNOWNS_AND_CROSS_BIASES = 18;
120
121 /**
122 * Number of unknowns for the general case when G-dependent cross
123 * biases are being estimated.
124 */
125 public static final int GENERAL_UNKNOWNS_AND_CROSS_BIASES = 21;
126
127 /**
128 * Number of unknowns when common z-axis is assumed for both
129 * the accelerometer and gyroscope when G-dependent cross biases
130 * are not being estimated.
131 */
132 public static final int COMMON_Z_AXIS_UNKNOWNS = 9;
133
134 /**
135 * Number of unknowns for the general case when G-dependent cross
136 * biases are not being estimated.
137 */
138 public static final int GENERAL_UNKNOWNS = 12;
139
140 /**
141 * Required minimum number of sequences when common z-axis is assumed
142 * and G-dependent cross biases are being estimated.
143 */
144 public static final int MINIMUM_SEQUENCES_COMMON_Z_AXIS_AND_CROSS_BIASES =
145 COMMON_Z_AXIS_UNKNOWNS_AND_CROSS_BIASES + 1;
146
147 /**
148 * Required minimum number of sequences for the general case and
149 * G-dependent cross biases are being estimated.
150 */
151 public static final int MINIMUM_SEQUENCES_GENERAL_AND_CROSS_BIASES =
152 GENERAL_UNKNOWNS_AND_CROSS_BIASES + 1;
153
154 /**
155 * Required minimum number of sequences when common z-axis is assumed
156 * and G-dependent cross biases are being ignored.
157 */
158 public static final int MINIMUM_SEQUENCES_COMMON_Z_AXIS = COMMON_Z_AXIS_UNKNOWNS + 1;
159
160 /**
161 * Required minimum number of sequences for the general case and
162 * G-dependent cross biases are being ignored.
163 */
164 public static final int MINIMUM_SEQUENCES_GENERAL = GENERAL_UNKNOWNS + 1;
165
166 /**
167 * Levenberg-Marquardt fitter to find a non-linear solution.
168 */
169 private final LevenbergMarquardtMultiDimensionFitter fitter = new LevenbergMarquardtMultiDimensionFitter();
170
171 /**
172 * Known x-coordinate of accelerometer bias to be used to fix measured
173 * specific force and find cross biases introduced by the accelerometer.
174 * This is expressed in meters per squared second (m/s^2).
175 */
176 private double accelerometerBiasX;
177
178 /**
179 * Known y-coordinate of accelerometer bias to be used to fix measured
180 * specific force and find cross biases introduced by the accelerometer.
181 * This is expressed in meters per squared second (m/s^2).
182 */
183 private double accelerometerBiasY;
184
185 /**
186 * Known z-coordinate of accelerometer bias to be used to fix measured
187 * specific force and find cross biases introduced by the accelerometer.
188 * This is expressed in meters per squared second (m/s^2).
189 */
190 private double accelerometerBiasZ;
191
192 /**
193 * Known accelerometer x scaling factor to be used to fix measured
194 * specific force and find cross biases introduced by the accelerometer.
195 */
196 private double accelerometerSx;
197
198 /**
199 * Known accelerometer y scaling factor to be used to fix measured
200 * specific force and find cross biases introduced by the accelerometer.
201 */
202 private double accelerometerSy;
203
204 /**
205 * Known accelerometer z scaling factor to be used to fix measured
206 * specific force and find cross biases introduced by the accelerometer.
207 */
208 private double accelerometerSz;
209
210 /**
211 * Known accelerometer x-y cross coupling error to be used to fix measured
212 * specific force and find cross biases introduced by the accelerometer.
213 */
214 private double accelerometerMxy;
215
216 /**
217 * Know accelerometer x-z cross coupling error to be used to fix measured
218 * specific force and find cross biases introduced by the accelerometer.
219 */
220 private double accelerometerMxz;
221
222 /**
223 * Known accelerometer y-x cross coupling error to be used to fix measured
224 * specific force and find cross biases introduced by the accelerometer.
225 */
226 private double accelerometerMyx;
227
228 /**
229 * Known accelerometer y-z cross coupling error to be used to fix measured
230 * specific force and find cross biases introduced by the accelerometer.
231 */
232 private double accelerometerMyz;
233
234 /**
235 * Known accelerometer z-x cross coupling error to be used to fix measured
236 * specific force and find cross biases introduced by the accelerometer.
237 */
238 private double accelerometerMzx;
239
240 /**
241 * Known accelerometer z-y cross coupling error to be used to fix measured
242 * specific force and find cross biases introduced by the accelerometer.
243 */
244 private double accelerometerMzy;
245
246 /**
247 * Initial x-coordinate of gyroscope bias to be used to find a solution.
248 * This is expressed in radians per second (rad/s).
249 */
250 private double initialBiasX;
251
252 /**
253 * Initial y-coordinate of gyroscope bias to be used to find a solution.
254 * This is expressed in radians per second (rad/s).
255 */
256 private double initialBiasY;
257
258 /**
259 * Initial z-coordinate of gyroscope bias to be used to find a solution.
260 * This is expressed in radians per second (rad/s).
261 */
262 private double initialBiasZ;
263
264 /**
265 * Initial gyroscope x scaling factor.
266 */
267 private double initialSx;
268
269 /**
270 * Initial gyroscope y scaling factor.
271 */
272 private double initialSy;
273
274 /**
275 * Initial gyroscope z scaling factor.
276 */
277 private double initialSz;
278
279 /**
280 * Initial gyroscope x-y cross coupling error.
281 */
282 private double initialMxy;
283
284 /**
285 * Initial gyroscope x-z cross coupling error.
286 */
287 private double initialMxz;
288
289 /**
290 * Initial gyroscope y-x cross coupling error.
291 */
292 private double initialMyx;
293
294 /**
295 * Initial gyroscope y-z cross coupling error.
296 */
297 private double initialMyz;
298
299 /**
300 * Initial gyroscope z-x cross coupling error.
301 */
302 private double initialMzx;
303
304 /**
305 * Initial gyroscope z-y cross coupling error.
306 */
307 private double initialMzy;
308
309 /**
310 * Initial G-dependent cross biases introduced on the gyroscope by the
311 * specific forces sensed by the accelerometer.
312 */
313 private Matrix initialGg;
314
315 /**
316 * Contains a collection of sequences of timestamped body kinematics
317 * measurements taken at a given position where the device moves freely
318 * with different orientations.
319 */
320 private List<BodyKinematicsSequence<StandardDeviationTimedBodyKinematics>> sequences;
321
322 /**
323 * This flag indicates whether z-axis is assumed to be common for accelerometer
324 * and gyroscope.
325 * When enabled, this eliminates 3 variables from Mg matrix.
326 */
327 private boolean commonAxisUsed = DEFAULT_USE_COMMON_Z_AXIS;
328
329 /**
330 * This flag indicates whether G-dependent cross biases are being
331 * estimated or not.
332 * When enabled, this adds 9 variables from Gg matrix.
333 */
334 private boolean estimateGDependentCrossBiases = DEFAULT_ESTIMATE_G_DEPENDENT_CROSS_BIASES;
335
336 /**
337 * Listener to handle events raised by this calibrator.
338 */
339 private EasyGyroscopeCalibratorListener listener;
340
341 /**
342 * Estimated angular rate biases for each IMU axis expressed in radians per
343 * second (rad/s).
344 */
345 private double[] estimatedBiases;
346
347 /**
348 * Estimated gyroscope scale factors and cross coupling errors.
349 * This is the product of matrix Tg containing cross coupling errors and Kg
350 * containing scaling factors.
351 * So that:
352 * <pre>
353 * Mg = [sx mxy mxz] = Tg*Kg
354 * [myx sy myz]
355 * [mzx mzy sz ]
356 * </pre>
357 * Where:
358 * <pre>
359 * Kg = [sx 0 0 ]
360 * [0 sy 0 ]
361 * [0 0 sz]
362 * </pre>
363 * and
364 * <pre>
365 * Tg = [1 -alphaXy alphaXz ]
366 * [alphaYx 1 -alphaYz]
367 * [-alphaZx alphaZy 1 ]
368 * </pre>
369 * Hence:
370 * <pre>
371 * Mg = [sx mxy mxz] = Tg*Kg = [sx -sy * alphaXy sz * alphaXz ]
372 * [myx sy myz] [sx * alphaYx sy -sz * alphaYz]
373 * [mzx mzy sz ] [-sx * alphaZx sy * alphaZy sz ]
374 * </pre>
375 * This instance allows any 3x3 matrix however, typically alphaYx, alphaZx and alphaZy
376 * are considered to be zero if the gyroscope z-axis is assumed to be the same
377 * as the body z-axis. When this is assumed, myx = mzx = mzy = 0 and the Mg matrix
378 * becomes upper diagonal:
379 * <pre>
380 * Mg = [sx mxy mxz]
381 * [0 sy myz]
382 * [0 0 sz ]
383 * </pre>
384 * Values of this matrix are unit-less.
385 */
386 private Matrix estimatedMg;
387
388 /**
389 * Estimated G-dependent cross biases introduced on the gyroscope by the
390 * specific forces sensed by the accelerometer.
391 * This instance allows any 3x3 matrix.
392 */
393 private Matrix estimatedGg;
394
395 /**
396 * Estimated covariance matrix for estimated parameters.
397 */
398 private Matrix estimatedCovariance;
399
400 /**
401 * Estimated chi square value.
402 */
403 private double estimatedChiSq;
404
405 /**
406 * Estimated degrees of freedom of chi square value. Degrees of freedom is equal to the number of sampled data
407 * minus the number of estimated parameters.
408 */
409 private int estimatedChiSqDegreesOfFreedom;
410
411 /**
412 * Estimated reduced chi square value. This is equal to estimated chi square value divided by its degrees of
413 * freedom. Ideally this value should be close to 1.0.
414 */
415 private double estimatedReducedChiSq;
416
417 /**
418 * Estimated mean square error respect to provided measurements.
419 */
420 private double estimatedMse;
421
422 /**
423 * Estimated probability of finding a smaller chi square value expressed as a value between 0.0 and 1.0. The smaller
424 * the found chi square value is, the better the fit of the estimated parameters to the actual parameter. Thus, the
425 * smaller the chance of finding a smaller chi square value, then the better the estimated fit is.
426 */
427 private double estimatedP;
428
429 /**
430 * Estimated measure of quality of estimated fit as a value between 0.0 and 1.0. The larger the quality value is,
431 * the better the fit that has been estimated.
432 */
433 private double estimatedQ;
434
435 /**
436 * Indicates whether calibrator is running.
437 */
438 private boolean running;
439
440 /**
441 * Acceleration fixer.
442 */
443 private final AccelerationFixer accelerationFixer = new AccelerationFixer();
444
445 /**
446 * Index of current point being evaluated.
447 */
448 private int i;
449
450 /**
451 * Holds integrated rotation of a sequence.
452 */
453 private final Quaternion q = new Quaternion();
454
455 /**
456 * Contains a copy of input sequences where fixed body kinematics will be updated.
457 */
458 private List<BodyKinematicsSequence<StandardDeviationTimedBodyKinematics>> fixedSequences;
459
460 /**
461 * Contains measured specific force on a sample within a sequence
462 * expressed as a 3x1 matrix.
463 */
464 private Matrix measuredSpecificForce;
465
466 /**
467 * Contains fixed specific force for a sample within a sequence
468 * using provided accelerometer parameters and expressed as a 3x1
469 * matrix.
470 */
471 private Matrix trueSpecificForce;
472
473 /**
474 * Contains measured angular rate on a sample within a sequence
475 * expressed as a 3x1 matrix.
476 */
477 private Matrix measuredAngularRate;
478
479 /**
480 * Contains fixed specific force for a sample within a sequence
481 * using current parameters being estimated.
482 */
483 private Matrix trueAngularRate;
484
485 /**
486 * Internally holds cross-coupling errors during calibration.
487 */
488 private Matrix m;
489
490 /**
491 * Internally holds inverse of cross-coupling errors during calibration.
492 */
493 private Matrix invM;
494
495 /**
496 * Internally holds biases during calibration.
497 */
498 private Matrix b;
499
500 /**
501 * Internally hold g-dependent cross biases during calibration.
502 */
503 private Matrix g;
504
505 /**
506 * Internally holds angular rate bias due to g-dependent cross biases
507 */
508 private Matrix tmp;
509
510 /**
511 * Holds normalized gravity versor at the start of a sequence.
512 * This is used to compute rotations of gravity versor during
513 * calibration.
514 */
515 private final InhomogeneousPoint3D startPoint = new InhomogeneousPoint3D();
516
517 /**
518 * Holds normalized gravity versor at the end of a sequence.
519 * This is used to compute rotations of gravity versor during
520 * calibration.
521 */
522 private final InhomogeneousPoint3D endPoint = new InhomogeneousPoint3D();
523
524 /**
525 * Contains expected gravity versor obtained from measurements fixed
526 * using known accelerometer parameters.
527 * This is reused for memory efficiency purposes during calibration.
528 */
529 private final InhomogeneousPoint3D expectedEndPoint = new InhomogeneousPoint3D();
530
531 /**
532 * Array containing normalized gravity versor before (former 3 values)
533 * and after (latter 3 values) a given sequence.
534 * This is used during calibration.
535 */
536 private double[] point;
537
538 /**
539 * Constructor.
540 */
541 public EasyGyroscopeCalibrator() {
542 try {
543 initialGg = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
544 } catch (final WrongSizeException ignore) {
545 // never happens
546 }
547 }
548
549 /**
550 * Constructor.
551 *
552 * @param sequences collection of sequences containing timestamped body
553 * kinematics measurements.
554 * @param initialBias initial gyroscope bias to be used to find a solution.
555 * This must be 3x1 and is expressed in radians per
556 * second (rad/s).
557 * @param initialMg initial gyroscope scale factors and cross coupling
558 * errors matrix. Must be 3x3.
559 * @param initialGg initial gyroscope G-dependent cross biases
560 * introduced on the gyroscope by the specific forces
561 * sensed by the accelerometer. Must be 3x3.
562 * @throws IllegalArgumentException if any of the provided values does
563 * not have proper size.
564 */
565 public EasyGyroscopeCalibrator(
566 final List<BodyKinematicsSequence<StandardDeviationTimedBodyKinematics>> sequences,
567 final Matrix initialBias, final Matrix initialMg, final Matrix initialGg) {
568 this();
569 this.sequences = sequences;
570 try {
571 setInitialBias(initialBias);
572 setInitialMg(initialMg);
573 setInitialGg(initialGg);
574 } catch (final LockedException ignore) {
575 // never happens
576 }
577 }
578
579 /**
580 * Constructor.
581 *
582 * @param sequences collection of sequences containing timestamped body
583 * kinematics measurements.
584 * @param initialBias initial gyroscope bias to be used to find a solution.
585 * This must be 3x1 and is expressed in radians per
586 * second (rad/s).
587 * @param initialMg initial gyroscope scale factors and cross coupling
588 * errors matrix. Must be 3x3.
589 * @param initialGg initial gyroscope G-dependent cross biases
590 * introduced on the gyroscope by the specific forces
591 * sensed by the accelerometer. Must be 3x3.
592 * @param listener listener to handle events raised by this
593 * calibrator.
594 * @throws IllegalArgumentException if any of the provided values does
595 * not have proper size.
596 */
597 public EasyGyroscopeCalibrator(
598 final List<BodyKinematicsSequence<StandardDeviationTimedBodyKinematics>> sequences,
599 final Matrix initialBias, final Matrix initialMg, final Matrix initialGg,
600 final EasyGyroscopeCalibratorListener listener) {
601 this(sequences, initialBias, initialMg, initialGg);
602 this.listener = listener;
603 }
604
605 /**
606 * Constructor.
607 *
608 * @param sequences collection of sequences containing timestamped body
609 * kinematics measurements.
610 * @param initialBias initial gyroscope bias to be used to find a
611 * solution. This must have length 3 and is expressed
612 * in radians per second (rad/s).
613 * @param initialMg initial gyroscope scale factors and cross coupling
614 * errors matrix. Must be 3x3.
615 * @param initialGg initial gyroscope G-dependent cross biases
616 * introduced on the gyroscope by the specific forces
617 * sensed by the accelerometer. Must be 3x3.
618 * @throws IllegalArgumentException if any of the provided values does
619 * not have proper size.
620 */
621 public EasyGyroscopeCalibrator(
622 final List<BodyKinematicsSequence<StandardDeviationTimedBodyKinematics>> sequences,
623 final double[] initialBias, final Matrix initialMg, final Matrix initialGg) {
624 this();
625 this.sequences = sequences;
626 try {
627 setInitialBias(initialBias);
628 setInitialMg(initialMg);
629 setInitialGg(initialGg);
630 } catch (final LockedException ignore) {
631 // never happens
632 }
633 }
634
635 /**
636 * Constructor.
637 *
638 * @param sequences collection of sequences containing timestamped body
639 * kinematics measurements.
640 * @param initialBias initial gyroscope bias to be used to find a
641 * solution. This must have length 3 and is expressed
642 * in radians per second (rad/s).
643 * @param initialMg initial gyroscope scale factors and cross coupling
644 * errors matrix. Must be 3x3.
645 * @param initialGg initial gyroscope G-dependent cross biases
646 * introduced on the gyroscope by the specific forces
647 * sensed by the accelerometer. Must be 3x3.
648 * @param listener listener to handle events raised by this
649 * calibrator.
650 * @throws IllegalArgumentException if any of the provided values does
651 * not have proper size.
652 */
653 public EasyGyroscopeCalibrator(
654 final List<BodyKinematicsSequence<StandardDeviationTimedBodyKinematics>> sequences,
655 final double[] initialBias, final Matrix initialMg, final Matrix initialGg,
656 final EasyGyroscopeCalibratorListener listener) {
657 this(sequences, initialBias, initialMg, initialGg);
658 this.listener = listener;
659 }
660
661 /**
662 * Constructor.
663 *
664 * @param sequences collection of sequences containing timestamped body
665 * kinematics measurements.
666 * @param initialBias initial gyroscope bias to be used to find a
667 * solution. This must have length 3 and is expressed
668 * in radians per second (rad/s).
669 * @param initialMg initial gyroscope scale factors and cross coupling
670 * errors matrix. Must be 3x3.
671 * @param initialGg initial gyroscope G-dependent cross biases
672 * introduced on the gyroscope by the specific forces
673 * sensed by the accelerometer. Must be 3x3.
674 * @param accelerometerBias known accelerometer bias. This must
675 * have length 3 and is expressed in
676 * meters per squared second
677 * (m/s^2).
678 * @param accelerometerMa known accelerometer scale factors and
679 * cross coupling matrix. Must be 3x3.
680 * @throws IllegalArgumentException if any of the provided values does
681 * not have proper size.
682 */
683 public EasyGyroscopeCalibrator(
684 final List<BodyKinematicsSequence<StandardDeviationTimedBodyKinematics>> sequences,
685 final double[] initialBias, final Matrix initialMg, final Matrix initialGg,
686 final double[] accelerometerBias, final Matrix accelerometerMa) {
687 this(sequences, initialBias, initialMg, initialGg);
688 try {
689 setAccelerometerBias(accelerometerBias);
690 setAccelerometerMa(accelerometerMa);
691 } catch (final LockedException ignore) {
692 // never happens
693 }
694 }
695
696 /**
697 * Constructor.
698 *
699 * @param sequences collection of sequences containing timestamped body
700 * kinematics measurements.
701 * @param initialBias initial gyroscope bias to be used to find a
702 * solution. This must have length 3 and is expressed
703 * in radians per second (rad/s).
704 * @param initialMg initial gyroscope scale factors and cross coupling
705 * errors matrix. Must be 3x3.
706 * @param initialGg initial gyroscope G-dependent cross biases
707 * introduced on the gyroscope by the specific forces
708 * sensed by the accelerometer. Must be 3x3.
709 * @param accelerometerBias known accelerometer bias. This must
710 * have length 3 and is expressed in
711 * meters per squared second
712 * (m/s^2).
713 * @param accelerometerMa known accelerometer scale factors and
714 * cross coupling matrix. Must be 3x3.
715 * @param listener listener to handle events raised by this
716 * calibrator.
717 * @throws IllegalArgumentException if any of the provided values does
718 * not have proper size.
719 */
720 public EasyGyroscopeCalibrator(
721 final List<BodyKinematicsSequence<StandardDeviationTimedBodyKinematics>> sequences,
722 final double[] initialBias, final Matrix initialMg, final Matrix initialGg,
723 final double[] accelerometerBias, final Matrix accelerometerMa,
724 final EasyGyroscopeCalibratorListener listener) {
725 this(sequences, initialBias, initialMg, initialGg, accelerometerBias, accelerometerMa);
726 this.listener = listener;
727 }
728
729 /**
730 * Constructor.
731 *
732 * @param sequences collection of sequences containing timestamped body
733 * kinematics measurements.
734 * @param initialBias initial gyroscope bias to be used to find a
735 * solution. This must be 3x1 and is expressed
736 * in radians per second (rad/s).
737 * @param initialMg initial gyroscope scale factors and cross coupling
738 * errors matrix. Must be 3x3.
739 * @param initialGg initial gyroscope G-dependent cross biases
740 * introduced on the gyroscope by the specific forces
741 * sensed by the accelerometer. Must be 3x3.
742 * @param accelerometerBias known accelerometer bias. This must be 3x1
743 * and is expressed in meters per squared
744 * second (m/s^2).
745 * @param accelerometerMa known accelerometer scale factors and
746 * cross coupling matrix. Must be 3x3.
747 * @throws IllegalArgumentException if any of the provided values does
748 * not have proper size.
749 */
750 public EasyGyroscopeCalibrator(
751 final List<BodyKinematicsSequence<StandardDeviationTimedBodyKinematics>> sequences,
752 final Matrix initialBias, final Matrix initialMg, final Matrix initialGg,
753 final Matrix accelerometerBias, final Matrix accelerometerMa) {
754 this(sequences, initialBias, initialMg, initialGg);
755 try {
756 setAccelerometerBias(accelerometerBias);
757 setAccelerometerMa(accelerometerMa);
758 } catch (final LockedException ignore) {
759 // never happens
760 }
761 }
762
763 /**
764 * Constructor.
765 *
766 * @param sequences collection of sequences containing timestamped body
767 * kinematics measurements.
768 * @param initialBias initial gyroscope bias to be used to find a
769 * solution. This must be 3x1 and is expressed
770 * in radians per second (rad/s).
771 * @param initialMg initial gyroscope scale factors and cross coupling
772 * errors matrix. Must be 3x3.
773 * @param initialGg initial gyroscope G-dependent cross biases
774 * introduced on the gyroscope by the specific forces
775 * sensed by the accelerometer. Must be 3x3.
776 * @param accelerometerBias known accelerometer bias. This must be 3x1
777 * and is expressed in meters per squared
778 * second (m/s^2).
779 * @param accelerometerMa known accelerometer scale factors and
780 * cross coupling matrix. Must be 3x3.
781 * @param listener listener to handle events raised by this
782 * calibrator.
783 * @throws IllegalArgumentException if any of the provided values does
784 * not have proper size.
785 */
786 public EasyGyroscopeCalibrator(
787 final List<BodyKinematicsSequence<StandardDeviationTimedBodyKinematics>> sequences,
788 final Matrix initialBias, final Matrix initialMg, final Matrix initialGg, final Matrix accelerometerBias,
789 final Matrix accelerometerMa, final EasyGyroscopeCalibratorListener listener) {
790 this(sequences, initialBias, initialMg, initialGg, accelerometerBias, accelerometerMa);
791 this.listener = listener;
792 }
793
794 /**
795 * Constructor.
796 *
797 * @param sequences collection of sequences containing timestamped body
798 * kinematics measurements.
799 * @param commonAxisUsed indicates whether z-axis is
800 * assumed to be common for
801 * accelerometer and gyroscope.
802 * @param estimateGDependentCrossBiases true if G-dependent cross biases
803 * will be estimated, false
804 * otherwise.
805 * @param initialBias initial gyroscope bias to be used to find a
806 * solution. This must be 3x1 and is expressed
807 * in radians per second (rad/s).
808 * @param initialMg initial gyroscope scale factors and cross coupling
809 * errors matrix. Must be 3x3.
810 * @param initialGg initial gyroscope G-dependent cross biases
811 * introduced on the gyroscope by the specific forces
812 * sensed by the accelerometer. Must be 3x3.
813 * @throws IllegalArgumentException if any of the provided values does
814 * not have proper size.
815 */
816 public EasyGyroscopeCalibrator(
817 final List<BodyKinematicsSequence<StandardDeviationTimedBodyKinematics>> sequences,
818 final boolean commonAxisUsed, final boolean estimateGDependentCrossBiases,
819 final Matrix initialBias, final Matrix initialMg, final Matrix initialGg) {
820 this(sequences, initialBias, initialMg, initialGg);
821 this.commonAxisUsed = commonAxisUsed;
822 this.estimateGDependentCrossBiases = estimateGDependentCrossBiases;
823 }
824
825 /**
826 * Constructor.
827 *
828 * @param sequences collection of sequences containing timestamped body
829 * kinematics measurements.
830 * @param commonAxisUsed indicates whether z-axis is
831 * assumed to be common for
832 * accelerometer and gyroscope.
833 * @param estimateGDependentCrossBiases true if G-dependent cross biases
834 * will be estimated, false
835 * otherwise.
836 * @param initialBias initial gyroscope bias to be used to find a
837 * solution. This must be 3x1 and is expressed
838 * in radians per second (rad/s).
839 * @param initialMg initial gyroscope scale factors and cross coupling
840 * errors matrix. Must be 3x3.
841 * @param initialGg initial gyroscope G-dependent cross biases
842 * introduced on the gyroscope by the specific forces
843 * sensed by the accelerometer. Must be 3x3.
844 * @param listener listener to handle events raised by this
845 * calibrator.
846 * @throws IllegalArgumentException if any of the provided values does
847 * not have proper size.
848 */
849 public EasyGyroscopeCalibrator(
850 final List<BodyKinematicsSequence<StandardDeviationTimedBodyKinematics>> sequences,
851 final boolean commonAxisUsed, final boolean estimateGDependentCrossBiases, final Matrix initialBias,
852 final Matrix initialMg, final Matrix initialGg, final EasyGyroscopeCalibratorListener listener) {
853 this(sequences, commonAxisUsed, estimateGDependentCrossBiases, initialBias, initialMg, initialGg);
854 this.listener = listener;
855 }
856
857 /**
858 * Constructor.
859 *
860 * @param sequences collection of sequences containing timestamped body
861 * kinematics measurements.
862 * @param commonAxisUsed indicates whether z-axis is
863 * assumed to be common for
864 * accelerometer and gyroscope.
865 * @param estimateGDependentCrossBiases true if G-dependent cross biases
866 * will be estimated, false
867 * otherwise.
868 * @param initialBias initial gyroscope bias to be used to find a
869 * solution. This must have length 3 and is expressed
870 * in radians per second (rad/s).
871 * @param initialMg initial gyroscope scale factors and cross coupling
872 * errors matrix. Must be 3x3.
873 * @param initialGg initial gyroscope G-dependent cross biases
874 * introduced on the gyroscope by the specific forces
875 * sensed by the accelerometer. Must be 3x3.
876 * @throws IllegalArgumentException if any of the provided values does
877 * not have proper size.
878 */
879 public EasyGyroscopeCalibrator(
880 final List<BodyKinematicsSequence<StandardDeviationTimedBodyKinematics>> sequences,
881 final boolean commonAxisUsed, final boolean estimateGDependentCrossBiases, final double[] initialBias,
882 final Matrix initialMg, final Matrix initialGg) {
883 this(sequences, initialBias, initialMg, initialGg);
884 this.commonAxisUsed = commonAxisUsed;
885 this.estimateGDependentCrossBiases = estimateGDependentCrossBiases;
886 }
887
888 /**
889 * Constructor.
890 *
891 * @param sequences collection of sequences containing timestamped body
892 * kinematics measurements.
893 * @param commonAxisUsed indicates whether z-axis is
894 * assumed to be common for
895 * accelerometer and gyroscope.
896 * @param estimateGDependentCrossBiases true if G-dependent cross biases
897 * will be estimated, false
898 * otherwise.
899 * @param initialBias initial gyroscope bias to be used to find a
900 * solution. This must have length 3 and is expressed
901 * in radians per second (rad/s).
902 * @param initialMg initial gyroscope scale factors and cross coupling
903 * errors matrix. Must be 3x3.
904 * @param initialGg initial gyroscope G-dependent cross biases
905 * introduced on the gyroscope by the specific forces
906 * sensed by the accelerometer. Must be 3x3.
907 * @param listener listener to handle events raised by this
908 * calibrator.
909 * @throws IllegalArgumentException if any of the provided values does
910 * not have proper size.
911 */
912 public EasyGyroscopeCalibrator(
913 final List<BodyKinematicsSequence<StandardDeviationTimedBodyKinematics>> sequences,
914 final boolean commonAxisUsed, final boolean estimateGDependentCrossBiases, final double[] initialBias,
915 final Matrix initialMg, final Matrix initialGg, final EasyGyroscopeCalibratorListener listener) {
916 this(sequences, commonAxisUsed, estimateGDependentCrossBiases, initialBias, initialMg, initialGg);
917 this.listener = listener;
918 }
919
920 /**
921 * Constructor.
922 *
923 * @param sequences collection of sequences containing timestamped body
924 * kinematics measurements.
925 * @param commonAxisUsed indicates whether z-axis is
926 * assumed to be common for
927 * accelerometer and gyroscope.
928 * @param estimateGDependentCrossBiases true if G-dependent cross biases
929 * will be estimated, false
930 * otherwise.
931 * @param initialBias initial gyroscope bias to be used to find a
932 * solution. This must have length 3 and is expressed
933 * in radians per second (rad/s).
934 * @param initialMg initial gyroscope scale factors and cross coupling
935 * errors matrix. Must be 3x3.
936 * @param initialGg initial gyroscope G-dependent cross biases
937 * introduced on the gyroscope by the specific forces
938 * sensed by the accelerometer. Must be 3x3.
939 * @param accelerometerBias known accelerometer bias. This
940 * must have length 3 and is
941 * expressed in meters per squared
942 * second (m/s^2).
943 * @param accelerometerMa known accelerometer scale factors
944 * and cross coupling matrix. Must
945 * be 3x3.
946 * @throws IllegalArgumentException if any of the provided values does
947 * not have proper size.
948 */
949 public EasyGyroscopeCalibrator(
950 final List<BodyKinematicsSequence<StandardDeviationTimedBodyKinematics>> sequences,
951 final boolean commonAxisUsed, final boolean estimateGDependentCrossBiases, final double[] initialBias,
952 final Matrix initialMg, final Matrix initialGg, final double[] accelerometerBias,
953 final Matrix accelerometerMa) {
954 this(sequences, initialBias, initialMg, initialGg, accelerometerBias, accelerometerMa);
955 this.commonAxisUsed = commonAxisUsed;
956 this.estimateGDependentCrossBiases = estimateGDependentCrossBiases;
957 }
958
959 /**
960 * Constructor.
961 *
962 * @param sequences collection of sequences containing timestamped body
963 * kinematics measurements.
964 * @param commonAxisUsed indicates whether z-axis is
965 * assumed to be common for
966 * accelerometer and gyroscope.
967 * @param estimateGDependentCrossBiases true if G-dependent cross biases
968 * will be estimated, false
969 * otherwise.
970 * @param initialBias initial gyroscope bias to be used to find a
971 * solution. This must have length 3 and is expressed
972 * in radians per second (rad/s).
973 * @param initialMg initial gyroscope scale factors and cross coupling
974 * errors matrix. Must be 3x3.
975 * @param initialGg initial gyroscope G-dependent cross biases
976 * introduced on the gyroscope by the specific forces
977 * sensed by the accelerometer. Must be 3x3.
978 * @param accelerometerBias known accelerometer bias. This
979 * must have length 3 and is
980 * expressed in meters per squared
981 * second (m/s^2).
982 * @param accelerometerMa known accelerometer scale factors
983 * and cross coupling matrix. Must
984 * be 3x3.
985 * @param listener listener to handle events raised by this
986 * calibrator.
987 * @throws IllegalArgumentException if any of the provided values does
988 * not have proper size.
989 */
990 public EasyGyroscopeCalibrator(
991 final List<BodyKinematicsSequence<StandardDeviationTimedBodyKinematics>> sequences,
992 final boolean commonAxisUsed, final boolean estimateGDependentCrossBiases, final double[] initialBias,
993 final Matrix initialMg, final Matrix initialGg, final double[] accelerometerBias,
994 final Matrix accelerometerMa, final EasyGyroscopeCalibratorListener listener) {
995 this(sequences, commonAxisUsed, estimateGDependentCrossBiases, initialBias, initialMg, initialGg,
996 accelerometerBias, accelerometerMa);
997 this.listener = listener;
998 }
999
1000 /**
1001 * Constructor.
1002 *
1003 * @param sequences collection of sequences containing timestamped body
1004 * kinematics measurements.
1005 * @param commonAxisUsed indicates whether z-axis is
1006 * assumed to be common for
1007 * accelerometer and gyroscope.
1008 * @param estimateGDependentCrossBiases true if G-dependent cross biases
1009 * will be estimated, false
1010 * otherwise.
1011 * @param initialBias initial gyroscope bias to be used to find a
1012 * solution. This must be 3x1 and is expressed
1013 * in radians per second (rad/s).
1014 * @param initialMg initial gyroscope scale factors and cross coupling
1015 * errors matrix. Must be 3x3.
1016 * @param initialGg initial gyroscope G-dependent cross biases
1017 * introduced on the gyroscope by the specific forces
1018 * sensed by the accelerometer. Must be 3x3.
1019 * @param accelerometerBias known accelerometer bias. This
1020 * must have length 3 and is
1021 * expressed in meters per squared
1022 * second (m/s^2).
1023 * @param accelerometerMa known accelerometer scale factors
1024 * and cross coupling matrix. Must
1025 * be 3x3.
1026 * @throws IllegalArgumentException if any of the provided values does
1027 * not have proper size.
1028 */
1029 public EasyGyroscopeCalibrator(
1030 final List<BodyKinematicsSequence<StandardDeviationTimedBodyKinematics>> sequences,
1031 final boolean commonAxisUsed, final boolean estimateGDependentCrossBiases, final Matrix initialBias,
1032 final Matrix initialMg, final Matrix initialGg, final Matrix accelerometerBias,
1033 final Matrix accelerometerMa) {
1034 this(sequences, initialBias, initialMg, initialGg, accelerometerBias, accelerometerMa);
1035 this.commonAxisUsed = commonAxisUsed;
1036 this.estimateGDependentCrossBiases = estimateGDependentCrossBiases;
1037 }
1038
1039 /**
1040 * Constructor.
1041 *
1042 * @param sequences collection of sequences containing timestamped body
1043 * kinematics measurements.
1044 * @param commonAxisUsed indicates whether z-axis is
1045 * assumed to be common for
1046 * accelerometer and gyroscope.
1047 * @param estimateGDependentCrossBiases true if G-dependent cross biases
1048 * will be estimated, false
1049 * otherwise.
1050 * @param initialBias initial gyroscope bias to be used to find a
1051 * solution. This must be 3x1 and is expressed
1052 * in radians per second (rad/s).
1053 * @param initialMg initial gyroscope scale factors and cross coupling
1054 * errors matrix. Must be 3x3.
1055 * @param initialGg initial gyroscope G-dependent cross biases
1056 * introduced on the gyroscope by the specific forces
1057 * sensed by the accelerometer. Must be 3x3.
1058 * @param accelerometerBias known accelerometer bias. This
1059 * must have length 3 and is
1060 * expressed in meters per squared
1061 * second (m/s^2).
1062 * @param accelerometerMa known accelerometer scale factors
1063 * and cross coupling matrix. Must
1064 * be 3x3.
1065 * @param listener listener to handle events raised by this
1066 * calibrator.
1067 * @throws IllegalArgumentException if any of the provided values does
1068 * not have proper size.
1069 */
1070 public EasyGyroscopeCalibrator(
1071 final List<BodyKinematicsSequence<StandardDeviationTimedBodyKinematics>> sequences,
1072 final boolean commonAxisUsed, final boolean estimateGDependentCrossBiases, final Matrix initialBias,
1073 final Matrix initialMg, final Matrix initialGg, final Matrix accelerometerBias,
1074 final Matrix accelerometerMa, final EasyGyroscopeCalibratorListener listener) {
1075 this(sequences, commonAxisUsed, estimateGDependentCrossBiases, initialBias, initialMg, initialGg,
1076 accelerometerBias, accelerometerMa);
1077 this.listener = listener;
1078 }
1079
1080 /**
1081 * Gets known x-coordinate of accelerometer bias to be used to fix
1082 * measured specific force and find cross biases introduced by the
1083 * accelerometer.
1084 * This is expressed in meters per squared second (m/s^2).
1085 *
1086 * @return known x-coordinate of accelerometer bias.
1087 */
1088 @Override
1089 public double getAccelerometerBiasX() {
1090 return accelerometerBiasX;
1091 }
1092
1093 /**
1094 * Sets known x-coordinate of accelerometer bias to be used to fix
1095 * measured specific force and find cross biases introduced by the
1096 * accelerometer.
1097 * This is expressed in meters per squared second (m/s^2).
1098 *
1099 * @param accelerometerBiasX known x-coordinate of accelerometer bias.
1100 * @throws LockedException if calibrator is currently running.
1101 */
1102 @Override
1103 public void setAccelerometerBiasX(final double accelerometerBiasX) throws LockedException {
1104 if (running) {
1105 throw new LockedException();
1106 }
1107 this.accelerometerBiasX = accelerometerBiasX;
1108 }
1109
1110 /**
1111 * Gets known y-coordinate of accelerometer bias to be used to fix
1112 * measured specific force and find cross biases introduced by the
1113 * accelerometer.
1114 * This is expressed in meters per squared second (m/s^2).
1115 *
1116 * @return known y-coordinate of accelerometer bias.
1117 */
1118 @Override
1119 public double getAccelerometerBiasY() {
1120 return accelerometerBiasY;
1121 }
1122
1123 /**
1124 * Sets known y-coordinate of accelerometer bias to be used to fix
1125 * measured specific force and find cross biases introduced by the
1126 * accelerometer.
1127 * This is expressed in meters per squared second (m/s^2).
1128 *
1129 * @param accelerometerBiasY known y-coordinate of accelerometer bias.
1130 * @throws LockedException if calibrator is currently running.
1131 */
1132 @Override
1133 public void setAccelerometerBiasY(final double accelerometerBiasY) throws LockedException {
1134 if (running) {
1135 throw new LockedException();
1136 }
1137 this.accelerometerBiasY = accelerometerBiasY;
1138 }
1139
1140 /**
1141 * Gets known z-coordinate of accelerometer bias to be used to fix
1142 * measured specific force and find cross biases introduced by the
1143 * accelerometer.
1144 * This is expressed in meters per squared second (m/s^2).
1145 *
1146 * @return known z-coordinate of accelerometer bias.
1147 */
1148 @Override
1149 public double getAccelerometerBiasZ() {
1150 return accelerometerBiasZ;
1151 }
1152
1153 /**
1154 * Sets known z-coordinate of accelerometer bias to be used to fix
1155 * measured specific force and find cross biases introduced by the
1156 * accelerometer.
1157 * This is expressed in meters per squared second (m/s^2).
1158 *
1159 * @param accelerometerBiasZ known z-coordinate of accelerometer bias.
1160 * @throws LockedException if calibrator is currently running.
1161 */
1162 @Override
1163 public void setAccelerometerBiasZ(final double accelerometerBiasZ) throws LockedException {
1164 if (running) {
1165 throw new LockedException();
1166 }
1167 this.accelerometerBiasZ = accelerometerBiasZ;
1168 }
1169
1170 /**
1171 * Gets known x-coordinate of accelerometer bias to be used to fix
1172 * measured specific force and find cross biases introduced by the
1173 * accelerometer.
1174 *
1175 * @return known x-coordinate of accelerometer bias.
1176 */
1177 @Override
1178 public Acceleration getAccelerometerBiasXAsAcceleration() {
1179 return new Acceleration(accelerometerBiasX, AccelerationUnit.METERS_PER_SQUARED_SECOND);
1180 }
1181
1182 /**
1183 * Gets known x-coordinate of accelerometer bias to be used to fix
1184 * measured specific force and find cross biases introduced by the
1185 * accelerometer.
1186 *
1187 * @param result instance where result data will be stored.
1188 */
1189 @Override
1190 public void getAccelerometerBiasXAsAcceleration(final Acceleration result) {
1191 result.setValue(accelerometerBiasX);
1192 result.setUnit(AccelerationUnit.METERS_PER_SQUARED_SECOND);
1193 }
1194
1195 /**
1196 * Sets known x-coordinate of accelerometer bias to be used to fix
1197 * measured specific force and find cross biases introduced by the
1198 * accelerometer.
1199 *
1200 * @param accelerometerBiasX x-coordinate of accelerometer bias.
1201 * @throws LockedException if calibrator is currently running.
1202 */
1203 @Override
1204 public void setAccelerometerBiasX(final Acceleration accelerometerBiasX) throws LockedException {
1205 if (running) {
1206 throw new LockedException();
1207 }
1208 this.accelerometerBiasX = convertAcceleration(accelerometerBiasX);
1209 }
1210
1211 /**
1212 * Gets known y-coordinate of accelerometer bias to be used to fix
1213 * measured specific force and find cross biases introduced by the
1214 * accelerometer.
1215 *
1216 * @return known y-coordinate of accelerometer bias.
1217 */
1218 @Override
1219 public Acceleration getAccelerometerBiasYAsAcceleration() {
1220 return new Acceleration(accelerometerBiasY, AccelerationUnit.METERS_PER_SQUARED_SECOND);
1221 }
1222
1223 /**
1224 * Gets known y-coordinate of accelerometer bias to be used to fix
1225 * measured specific force and find cross biases introduced by the
1226 * accelerometer.
1227 *
1228 * @param result instance where result data will be stored.
1229 */
1230 @Override
1231 public void getAccelerometerBiasYAsAcceleration(final Acceleration result) {
1232 result.setValue(accelerometerBiasY);
1233 result.setUnit(AccelerationUnit.METERS_PER_SQUARED_SECOND);
1234 }
1235
1236 /**
1237 * Sets known y-coordinate of accelerometer bias to be used to fix
1238 * measured specific force and find cross biases introduced by the
1239 * accelerometer.
1240 *
1241 * @param accelerometerBiasY y-coordinate of accelerometer bias.
1242 * @throws LockedException if calibrator is currently running.
1243 */
1244 @Override
1245 public void setAccelerometerBiasY(final Acceleration accelerometerBiasY) throws LockedException {
1246 if (running) {
1247 throw new LockedException();
1248 }
1249 this.accelerometerBiasY = convertAcceleration(accelerometerBiasY);
1250 }
1251
1252 /**
1253 * Gets known z-coordinate of accelerometer bias to be used to fix
1254 * measured specific force and find cross biases introduced by the
1255 * accelerometer.
1256 *
1257 * @return known z-coordinate of accelerometer bias.
1258 */
1259 @Override
1260 public Acceleration getAccelerometerBiasZAsAcceleration() {
1261 return new Acceleration(accelerometerBiasZ, AccelerationUnit.METERS_PER_SQUARED_SECOND);
1262 }
1263
1264 /**
1265 * Gets known z-coordinate of accelerometer bias to be used to fix
1266 * measured specific force and find cross biases introduced by the
1267 * accelerometer.
1268 *
1269 * @param result instance where result data will be stored.
1270 */
1271 @Override
1272 public void getAccelerometerBiasZAsAcceleration(final Acceleration result) {
1273 result.setValue(accelerometerBiasZ);
1274 result.setUnit(AccelerationUnit.METERS_PER_SQUARED_SECOND);
1275 }
1276
1277 /**
1278 * Sets known z-coordinate of accelerometer bias to be used to fix
1279 * measured specific force and find cross biases introduced by the
1280 * accelerometer.
1281 *
1282 * @param accelerometerBiasZ z-coordinate of accelerometer bias.
1283 * @throws LockedException if calibrator is currently running.
1284 */
1285 @Override
1286 public void setAccelerometerBiasZ(final Acceleration accelerometerBiasZ) throws LockedException {
1287 if (running) {
1288 throw new LockedException();
1289 }
1290 this.accelerometerBiasZ = convertAcceleration(accelerometerBiasZ);
1291 }
1292
1293 /**
1294 * Sets known accelerometer bias to be used to fix measured specific
1295 * force and find cross biases introduced by the accelerometer.
1296 * This is expressed in meters per squared second (m/s^2).
1297 *
1298 * @param accelerometerBiasX x-coordinate of accelerometer bias.
1299 * @param accelerometerBiasY y-coordinate of accelerometer bias.
1300 * @param accelerometerBiasZ z-coordinate of accelerometer bias.
1301 * @throws LockedException if calibrator is currently running.
1302 */
1303 @Override
1304 public void setAccelerometerBias(
1305 final double accelerometerBiasX, final double accelerometerBiasY, final double accelerometerBiasZ)
1306 throws LockedException {
1307 if (running) {
1308 throw new LockedException();
1309 }
1310
1311 this.accelerometerBiasX = accelerometerBiasX;
1312 this.accelerometerBiasY = accelerometerBiasY;
1313 this.accelerometerBiasZ = accelerometerBiasZ;
1314 }
1315
1316 /**
1317 * Sets known accelerometer bias to be used to fix measured specific
1318 * force and find cross biases introduced by the accelerometer.
1319 *
1320 * @param accelerometerBiasX x-coordinate of accelerometer bias.
1321 * @param accelerometerBiasY y-coordinate of accelerometer bias.
1322 * @param accelerometerBiasZ z-coordinate of accelerometer bias.
1323 * @throws LockedException if calibrator is currently running.
1324 */
1325 @Override
1326 public void setAccelerometerBias(
1327 final Acceleration accelerometerBiasX,
1328 final Acceleration accelerometerBiasY,
1329 final Acceleration accelerometerBiasZ) throws LockedException {
1330 if (running) {
1331 throw new LockedException();
1332 }
1333
1334 this.accelerometerBiasX = convertAcceleration(accelerometerBiasX);
1335 this.accelerometerBiasY = convertAcceleration(accelerometerBiasY);
1336 this.accelerometerBiasZ = convertAcceleration(accelerometerBiasZ);
1337 }
1338
1339 /**
1340 * Gets known accelerometer bias to be used to fix measured specific
1341 * force and find cross biases introduced by the accelerometer.
1342 * This is expressed in meters per squared second (m/s^2).
1343 *
1344 * @return known accelerometer bias.
1345 */
1346 @Override
1347 public double[] getAccelerometerBias() {
1348 final var result = new double[BodyKinematics.COMPONENTS];
1349 getAccelerometerBias(result);
1350 return result;
1351 }
1352
1353 /**
1354 * Gets known accelerometer bias to be used to fix measured specific
1355 * force and find cross biases introduced by the accelerometer.
1356 * This is expressed in meters per squared second (m/s^2).
1357 *
1358 * @param result instance where result data will be copied to.
1359 * @throws IllegalArgumentException if provided array does not have
1360 * length 3.
1361 */
1362 @Override
1363 public void getAccelerometerBias(final double[] result) {
1364 if (result.length != BodyKinematics.COMPONENTS) {
1365 throw new IllegalArgumentException();
1366 }
1367
1368 result[0] = accelerometerBiasX;
1369 result[1] = accelerometerBiasY;
1370 result[2] = accelerometerBiasZ;
1371 }
1372
1373 /**
1374 * Sets known accelerometer bias to be used to fix measured specific
1375 * force and find cross biases introduced by the accelerometer.
1376 * This is expressed in meters per squared second (m/s^2).
1377 *
1378 * @param accelerometerBias known accelerometer bias.
1379 * @throws LockedException if calibrator is currently running.
1380 * @throws IllegalArgumentException if provided array does not have
1381 * length 3.
1382 */
1383 @Override
1384 public void setAccelerometerBias(final double[] accelerometerBias) throws LockedException {
1385 if (running) {
1386 throw new LockedException();
1387 }
1388
1389 if (accelerometerBias.length != BodyKinematics.COMPONENTS) {
1390 throw new IllegalArgumentException();
1391 }
1392
1393 accelerometerBiasX = accelerometerBias[0];
1394 accelerometerBiasY = accelerometerBias[1];
1395 accelerometerBiasZ = accelerometerBias[2];
1396 }
1397
1398 /**
1399 * Gets known accelerometer bias to be used to fix measured specific
1400 * force and find cross biases introduced by the accelerometer.
1401 * This is expressed in meters per squared second (m/s^2).
1402 *
1403 * @return known accelerometer bias.
1404 */
1405 @Override
1406 public Matrix getAccelerometerBiasAsMatrix() {
1407 Matrix result;
1408 try {
1409 result = new Matrix(BodyKinematics.COMPONENTS, 1);
1410 getAccelerometerBiasAsMatrix(result);
1411 } catch (final WrongSizeException ignore) {
1412 // never happens
1413 result = null;
1414 }
1415 return result;
1416 }
1417
1418 /**
1419 * Gets known accelerometer bias to be used to fix measured specific
1420 * force and find cross biases introduced by the accelerometer.
1421 * This is expressed in meters per squared second (m/s^2).
1422 *
1423 * @param result instance where result data will be copied to.
1424 * @throws IllegalArgumentException if provided matrix is not 3x1.
1425 */
1426 @Override
1427 public void getAccelerometerBiasAsMatrix(final Matrix result) {
1428 if (result.getRows() != BodyKinematics.COMPONENTS || result.getColumns() != 1) {
1429 throw new IllegalArgumentException();
1430 }
1431 result.setElementAtIndex(0, accelerometerBiasX);
1432 result.setElementAtIndex(1, accelerometerBiasY);
1433 result.setElementAtIndex(2, accelerometerBiasZ);
1434 }
1435
1436 /**
1437 * Sets known accelerometer bias to be used to fix measured specific
1438 * force and find cross biases introduced by the accelerometer.
1439 * This is expressed in meters per squared second (m/s^2).
1440 *
1441 * @param accelerometerBias known accelerometer bias. Must be 3x1.
1442 * @throws LockedException if calibrator is currently running.
1443 * @throws IllegalArgumentException if provided matrix is not 3x1.
1444 */
1445 @Override
1446 public void setAccelerometerBias(final Matrix accelerometerBias) throws LockedException {
1447 if (running) {
1448 throw new LockedException();
1449 }
1450 if (accelerometerBias.getRows() != BodyKinematics.COMPONENTS || accelerometerBias.getColumns() != 1) {
1451 throw new IllegalArgumentException();
1452 }
1453
1454 accelerometerBiasX = accelerometerBias.getElementAtIndex(0);
1455 accelerometerBiasY = accelerometerBias.getElementAtIndex(1);
1456 accelerometerBiasZ = accelerometerBias.getElementAtIndex(2);
1457 }
1458
1459 /**
1460 * Gets known accelerometer x scaling factor to be used to fix measured
1461 * specific force and find cross biases introduced by the accelerometer.
1462 *
1463 * @return known accelerometer x scaling factor.
1464 */
1465 @Override
1466 public double getAccelerometerSx() {
1467 return accelerometerSx;
1468 }
1469
1470 /**
1471 * Sets known accelerometer x scaling factor to be used to fix measured
1472 * specific force and find cross biases introduced by the accelerometer.
1473 *
1474 * @param accelerometerSx known accelerometer x scaling factor.
1475 * @throws LockedException if calibrator is currently running.
1476 */
1477 @Override
1478 public void setAccelerometerSx(final double accelerometerSx) throws LockedException {
1479 if (running) {
1480 throw new LockedException();
1481 }
1482 this.accelerometerSx = accelerometerSx;
1483 }
1484
1485 /**
1486 * Gets known accelerometer y scaling factor to be used to fix measured
1487 * specific force and find cross biases introduced by the accelerometer.
1488 *
1489 * @return known accelerometer y scaling factor.
1490 */
1491 @Override
1492 public double getAccelerometerSy() {
1493 return accelerometerSy;
1494 }
1495
1496 /**
1497 * Sets known accelerometer y scaling factor to be used to fix measured
1498 * specific force and find cross biases introduced by the accelerometer.
1499 *
1500 * @param accelerometerSy known accelerometer y scaling factor.
1501 * @throws LockedException if calibrator is currently running.
1502 */
1503 @Override
1504 public void setAccelerometerSy(final double accelerometerSy) throws LockedException {
1505 if (running) {
1506 throw new LockedException();
1507 }
1508 this.accelerometerSy = accelerometerSy;
1509 }
1510
1511 /**
1512 * Gets known accelerometer z scaling factor to be used to fix measured
1513 * specific force and find cross biases introduced by the accelerometer.
1514 *
1515 * @return known accelerometer z scaling factor.
1516 */
1517 @Override
1518 public double getAccelerometerSz() {
1519 return accelerometerSz;
1520 }
1521
1522 /**
1523 * Sets known accelerometer z scaling factor to be used to fix measured
1524 * specific force and find cross biases introduced by the accelerometer.
1525 *
1526 * @param accelerometerSz known accelerometer z scaling factor.
1527 * @throws LockedException if calibrator is currently running.
1528 */
1529 @Override
1530 public void setAccelerometerSz(final double accelerometerSz) throws LockedException {
1531 if (running) {
1532 throw new LockedException();
1533 }
1534 this.accelerometerSz = accelerometerSz;
1535 }
1536
1537 /**
1538 * Gets known accelerometer x-y cross coupling error to be used to fix
1539 * measured specific force and find cross biases introduced by the
1540 * accelerometer.
1541 *
1542 * @return known accelerometer x-y cross coupling error.
1543 */
1544 @Override
1545 public double getAccelerometerMxy() {
1546 return accelerometerMxy;
1547 }
1548
1549 /**
1550 * Sets known accelerometer x-y cross coupling error to be used to fix
1551 * measured specific force and find cross biases introduced by the
1552 * accelerometer.
1553 *
1554 * @param accelerometerMxy known accelerometer x-y cross coupling error.
1555 * @throws LockedException if calibrator is currently running.
1556 */
1557 @Override
1558 public void setAccelerometerMxy(final double accelerometerMxy) throws LockedException {
1559 if (running) {
1560 throw new LockedException();
1561 }
1562 this.accelerometerMxy = accelerometerMxy;
1563 }
1564
1565 /**
1566 * Gets known accelerometer x-z cross coupling error to be used to fix
1567 * measured specific force and find cross biases introduced by the
1568 * accelerometer.
1569 *
1570 * @return known accelerometer x-z cross coupling error.
1571 */
1572 @Override
1573 public double getAccelerometerMxz() {
1574 return accelerometerMxz;
1575 }
1576
1577 /**
1578 * Sets known accelerometer x-z cross coupling error to be used to fix
1579 * measured specific force and find cross biases introduced by the
1580 * accelerometer.
1581 *
1582 * @param accelerometerMxz known accelerometer x-z cross coupling error.
1583 * @throws LockedException if calibrator is currently running.
1584 */
1585 @Override
1586 public void setAccelerometerMxz(final double accelerometerMxz) throws LockedException {
1587 if (running) {
1588 throw new LockedException();
1589 }
1590 this.accelerometerMxz = accelerometerMxz;
1591 }
1592
1593 /**
1594 * Gets known accelerometer y-x cross coupling error to be used to fix
1595 * measured specific force and find cross biases introduced by the
1596 * accelerometer.
1597 *
1598 * @return known accelerometer y-x cross coupling error.
1599 */
1600 @Override
1601 public double getAccelerometerMyx() {
1602 return accelerometerMyx;
1603 }
1604
1605 /**
1606 * Sets known accelerometer y-x cross coupling error to be used to fix
1607 * measured specific force and find cross biases introduced by the
1608 * accelerometer.
1609 *
1610 * @param accelerometerMyx known accelerometer y-x cross coupling
1611 * error.
1612 * @throws LockedException if calibrator is currently running.
1613 */
1614 @Override
1615 public void setAccelerometerMyx(final double accelerometerMyx) throws LockedException {
1616 if (running) {
1617 throw new LockedException();
1618 }
1619 this.accelerometerMyx = accelerometerMyx;
1620 }
1621
1622 /**
1623 * Gets known accelerometer y-z cross coupling error to be used to fix
1624 * measured specific force and find cross biases introduced by the
1625 * accelerometer.
1626 *
1627 * @return known accelerometer y-z cross coupling error.
1628 */
1629 @Override
1630 public double getAccelerometerMyz() {
1631 return accelerometerMyz;
1632 }
1633
1634 /**
1635 * Sets known accelerometer y-z cross coupling error to be used to fix
1636 * measured specific force and find cross biases introduced by the
1637 * accelerometer.
1638 *
1639 * @param accelerometerMyz known accelerometer y-z cross coupling
1640 * error.
1641 * @throws LockedException if calibrator is currently running.
1642 */
1643 @Override
1644 public void setAccelerometerMyz(final double accelerometerMyz) throws LockedException {
1645 if (running) {
1646 throw new LockedException();
1647 }
1648 this.accelerometerMyz = accelerometerMyz;
1649 }
1650
1651 /**
1652 * Gets known accelerometer z-x cross coupling error to be used to fix
1653 * measured specific force and find cross biases introduced by the
1654 * accelerometer.
1655 *
1656 * @return known accelerometer z-x cross coupling error.
1657 */
1658 @Override
1659 public double getAccelerometerMzx() {
1660 return accelerometerMzx;
1661 }
1662
1663 /**
1664 * Sets known accelerometer z-x cross coupling error to be used to fix
1665 * measured specific force and find cross biases introduced by the
1666 * accelerometer.
1667 *
1668 * @param accelerometerMzx known accelerometer z-x cross coupling
1669 * error.
1670 * @throws LockedException if calibrator is currently running.
1671 */
1672 @Override
1673 public void setAccelerometerMzx(final double accelerometerMzx) throws LockedException {
1674 if (running) {
1675 throw new LockedException();
1676 }
1677 this.accelerometerMzx = accelerometerMzx;
1678 }
1679
1680 /**
1681 * Gets known accelerometer z-y cross coupling error to be used to fix
1682 * measured specific force and find cross biases introduced by the
1683 * accelerometer.
1684 *
1685 * @return known accelerometer z-y cross coupling error.
1686 */
1687 @Override
1688 public double getAccelerometerMzy() {
1689 return accelerometerMzy;
1690 }
1691
1692 /**
1693 * Sets known accelerometer z-y cross coupling error to be used to fix
1694 * measured specific force and find cross biases introduced by the
1695 * accelerometer.
1696 *
1697 * @param accelerometerMzy known accelerometer z-y cross coupling
1698 * error.
1699 * @throws LockedException if calibrator is currently running.
1700 */
1701 @Override
1702 public void setAccelerometerMzy(final double accelerometerMzy) throws LockedException {
1703 if (running) {
1704 throw new LockedException();
1705 }
1706 this.accelerometerMzy = accelerometerMzy;
1707 }
1708
1709 /**
1710 * Sets known accelerometer scaling factors to be used to fix measured
1711 * specific force and find cross biases introduced by the
1712 * accelerometer.
1713 *
1714 * @param accelerometerSx known accelerometer x scaling factor.
1715 * @param accelerometerSy known accelerometer y scaling factor.
1716 * @param accelerometerSz known accelerometer z scaling factor.
1717 * @throws LockedException if calibrator is currently running.
1718 */
1719 @Override
1720 public void setAccelerometerScalingFactors(
1721 final double accelerometerSx, final double accelerometerSy, final double accelerometerSz)
1722 throws LockedException {
1723 if (running) {
1724 throw new LockedException();
1725 }
1726 this.accelerometerSx = accelerometerSx;
1727 this.accelerometerSy = accelerometerSy;
1728 this.accelerometerSz = accelerometerSz;
1729 }
1730
1731 /**
1732 * Sets known accelerometer cross coupling errors to be used to fix
1733 * measured specific force and find cross biases introduced by the
1734 * accelerometer.
1735 *
1736 * @param accelerometerMxy known accelerometer x-y cross coupling
1737 * error.
1738 * @param accelerometerMxz known accelerometer x-z cross coupling
1739 * error.
1740 * @param accelerometerMyx known accelerometer y-x cross coupling
1741 * error.
1742 * @param accelerometerMyz known accelerometer y-z cross coupling
1743 * error.
1744 * @param accelerometerMzx known accelerometer z-x cross coupling
1745 * error.
1746 * @param accelerometerMzy known accelerometer z-y cross coupling
1747 * error.
1748 * @throws LockedException if calibrator is currently running.
1749 */
1750 @Override
1751 public void setAccelerometerCrossCouplingErrors(
1752 final double accelerometerMxy, final double accelerometerMxz,
1753 final double accelerometerMyx, final double accelerometerMyz,
1754 final double accelerometerMzx, final double accelerometerMzy) throws LockedException {
1755 if (running) {
1756 throw new LockedException();
1757 }
1758 this.accelerometerMxy = accelerometerMxy;
1759 this.accelerometerMxz = accelerometerMxz;
1760 this.accelerometerMyx = accelerometerMyx;
1761 this.accelerometerMyz = accelerometerMyz;
1762 this.accelerometerMzx = accelerometerMzx;
1763 this.accelerometerMzy = accelerometerMzy;
1764 }
1765
1766 /**
1767 * Sets known accelerometer scaling factors and cross coupling errors
1768 * to be used to fix measured specific force and find cross biases
1769 * introduced by the accelerometer.
1770 *
1771 * @param accelerometerSx known accelerometer x scaling factor.
1772 * @param accelerometerSy known accelerometer y scaling factor.
1773 * @param accelerometerSz known accelerometer z scaling factor.
1774 * @param accelerometerMxy known accelerometer x-y cross coupling
1775 * error.
1776 * @param accelerometerMxz known accelerometer x-z cross coupling
1777 * error.
1778 * @param accelerometerMyx known accelerometer y-x cross coupling
1779 * error.
1780 * @param accelerometerMyz known accelerometer y-z cross coupling
1781 * error.
1782 * @param accelerometerMzx known accelerometer z-x cross coupling
1783 * error.
1784 * @param accelerometerMzy known accelerometer z-y cross coupling
1785 * error.
1786 * @throws LockedException if calibrator is currently running.
1787 */
1788 @Override
1789 public void setAccelerometerScalingFactorsAndCrossCouplingErrors(
1790 final double accelerometerSx, final double accelerometerSy,
1791 final double accelerometerSz, final double accelerometerMxy,
1792 final double accelerometerMxz, final double accelerometerMyx,
1793 final double accelerometerMyz, final double accelerometerMzx,
1794 final double accelerometerMzy) throws LockedException {
1795 if (running) {
1796 throw new LockedException();
1797 }
1798 setAccelerometerScalingFactors(accelerometerSx, accelerometerSy, accelerometerSz);
1799 setAccelerometerCrossCouplingErrors(accelerometerMxy, accelerometerMxz, accelerometerMyx,
1800 accelerometerMyz, accelerometerMzx, accelerometerMzy);
1801 }
1802
1803 /**
1804 * Gets known accelerometer scale factors and cross coupling
1805 * errors matrix.
1806 *
1807 * @return known accelerometer scale factors and cross coupling
1808 * errors matrix.
1809 */
1810 @Override
1811 public Matrix getAccelerometerMa() {
1812 Matrix result;
1813 try {
1814 result = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
1815 getAccelerometerMa(result);
1816 } catch (final WrongSizeException ignore) {
1817 // never happens
1818 result = null;
1819 }
1820 return result;
1821 }
1822
1823 /**
1824 * Gets known accelerometer scale factors and cross coupling
1825 * errors matrix.
1826 *
1827 * @param result instance where data will be stored.
1828 * @throws IllegalArgumentException if provided matrix is not 3x3.
1829 */
1830 @Override
1831 public void getAccelerometerMa(final Matrix result) {
1832 if (result.getRows() != BodyKinematics.COMPONENTS || result.getColumns() != BodyKinematics.COMPONENTS) {
1833 throw new IllegalArgumentException();
1834 }
1835 result.setElementAtIndex(0, accelerometerSx);
1836 result.setElementAtIndex(1, accelerometerMyx);
1837 result.setElementAtIndex(2, accelerometerMzx);
1838
1839 result.setElementAtIndex(3, accelerometerMxy);
1840 result.setElementAtIndex(4, accelerometerSy);
1841 result.setElementAtIndex(5, accelerometerMzy);
1842
1843 result.setElementAtIndex(6, accelerometerMxz);
1844 result.setElementAtIndex(7, accelerometerMyz);
1845 result.setElementAtIndex(8, accelerometerSz);
1846 }
1847
1848 /**
1849 * Sets known accelerometer scale factors and cross coupling
1850 * errors matrix.
1851 *
1852 * @param accelerometerMa known accelerometer scale factors and
1853 * cross coupling errors matrix. Must be 3x3.
1854 * @throws LockedException if calibrator is currently running.
1855 * @throws IllegalArgumentException if provided matrix is not 3x3.
1856 */
1857 @Override
1858 public void setAccelerometerMa(final Matrix accelerometerMa) throws LockedException {
1859 if (running) {
1860 throw new LockedException();
1861 }
1862 if (accelerometerMa.getRows() != BodyKinematics.COMPONENTS
1863 || accelerometerMa.getColumns() != BodyKinematics.COMPONENTS) {
1864 throw new IllegalArgumentException();
1865 }
1866
1867 accelerometerSx = accelerometerMa.getElementAtIndex(0);
1868 accelerometerMyx = accelerometerMa.getElementAtIndex(1);
1869 accelerometerMzx = accelerometerMa.getElementAtIndex(2);
1870
1871 accelerometerMxy = accelerometerMa.getElementAtIndex(3);
1872 accelerometerSy = accelerometerMa.getElementAtIndex(4);
1873 accelerometerMzy = accelerometerMa.getElementAtIndex(5);
1874
1875 accelerometerMxz = accelerometerMa.getElementAtIndex(6);
1876 accelerometerMyz = accelerometerMa.getElementAtIndex(7);
1877 accelerometerSz = accelerometerMa.getElementAtIndex(8);
1878 }
1879
1880 /**
1881 * Gets initial x-coordinate of gyroscope bias to be used to find
1882 * a solution.
1883 * This is expressed in radians per second (rad/s).
1884 *
1885 * @return initial x-coordinate of gyroscope bias.
1886 */
1887 public double getInitialBiasX() {
1888 return initialBiasX;
1889 }
1890
1891 /**
1892 * Sets initial x-coordinate of gyroscope bias to be used to find
1893 * a solution.
1894 * This is expressed in radians per second (rad/s).
1895 *
1896 * @param initialBiasX initial x-coordinate of gyroscope bias.
1897 * @throws LockedException if calibrator is currently running.
1898 */
1899 public void setInitialBiasX(final double initialBiasX) throws LockedException {
1900 if (running) {
1901 throw new LockedException();
1902 }
1903 this.initialBiasX = initialBiasX;
1904 }
1905
1906 /**
1907 * Gets initial y-coordinate of gyroscope bias to be used to find
1908 * a solution.
1909 * This is expressed in radians per second (rad/s).
1910 *
1911 * @return initial y-coordinate of gyroscope bias.
1912 */
1913 public double getInitialBiasY() {
1914 return initialBiasY;
1915 }
1916
1917 /**
1918 * Sets initial y-coordinate of gyroscope bias to be used to find
1919 * a solution.
1920 * This is expressed in radians per second (rad/s).
1921 *
1922 * @param initialBiasY initial y-coordinate of gyroscope bias.
1923 * @throws LockedException if calibrator is currently running.
1924 */
1925 public void setInitialBiasY(final double initialBiasY) throws LockedException {
1926 if (running) {
1927 throw new LockedException();
1928 }
1929 this.initialBiasY = initialBiasY;
1930 }
1931
1932 /**
1933 * Gets initial z-coordinate of gyroscope bias ot be used to find
1934 * a solution.
1935 * This is expressed in radians per second (rad/s).
1936 *
1937 * @return initial z-coordinate of gyroscope bias.
1938 */
1939 public double getInitialBiasZ() {
1940 return initialBiasZ;
1941 }
1942
1943 /**
1944 * Sets initial z-coordinate of gyroscope bias to be used to find
1945 * a solution.
1946 * This is expressed in radians per second (rad/s).
1947 *
1948 * @param initialBiasZ initial z-coordinate of gyroscope bias.
1949 * @throws LockedException if calibrator is currently running.
1950 */
1951 public void setInitialBiasZ(final double initialBiasZ) throws LockedException {
1952 if (running) {
1953 throw new LockedException();
1954 }
1955 this.initialBiasZ = initialBiasZ;
1956 }
1957
1958 /**
1959 * Gets initial x-coordinate of gyroscope bias to be used to find a
1960 * solution.
1961 *
1962 * @return initial x-coordinate of gyroscope bias.
1963 */
1964 public AngularSpeed getInitialBiasAngularSpeedX() {
1965 return new AngularSpeed(initialBiasX, AngularSpeedUnit.RADIANS_PER_SECOND);
1966 }
1967
1968 /**
1969 * Gets initial x-coordinate of gyroscope bias to be used to find a
1970 * solution.
1971 *
1972 * @param result instance where result data will be stored.
1973 */
1974 public void getInitialBiasAngularSpeedX(final AngularSpeed result) {
1975 result.setValue(initialBiasX);
1976 result.setUnit(AngularSpeedUnit.RADIANS_PER_SECOND);
1977 }
1978
1979 /**
1980 * Sets initial x-coordinate of gyroscope bias to be used to find a
1981 * solution.
1982 *
1983 * @param initialBiasX initial x-coordinate of gyroscope bias.
1984 * @throws LockedException if calibrator is currently running.
1985 */
1986 public void setInitialBiasX(final AngularSpeed initialBiasX) throws LockedException {
1987 if (running) {
1988 throw new LockedException();
1989 }
1990 this.initialBiasX = convertAngularSpeed(initialBiasX);
1991 }
1992
1993 /**
1994 * Gets initial y-coordinate of gyroscope bias to be used to find a
1995 * solution.
1996 *
1997 * @return initial y-coordinate of gyroscope bias.
1998 */
1999 public AngularSpeed getInitialBiasAngularSpeedY() {
2000 return new AngularSpeed(initialBiasY, AngularSpeedUnit.RADIANS_PER_SECOND);
2001 }
2002
2003 /**
2004 * Gets initial y-coordinate of gyroscope bias to be used to find a
2005 * solution.
2006 *
2007 * @param result instance where result data will be stored.
2008 */
2009 public void getInitialBiasAngularSpeedY(final AngularSpeed result) {
2010 result.setValue(initialBiasY);
2011 result.setUnit(AngularSpeedUnit.RADIANS_PER_SECOND);
2012 }
2013
2014 /**
2015 * Sets initial y-coordinate of gyroscope bias to be used to find a
2016 * solution.
2017 *
2018 * @param initialBiasY initial y-coordinate of gyroscope bias.
2019 * @throws LockedException if calibrator is currently running.
2020 */
2021 public void setInitialBiasY(final AngularSpeed initialBiasY) throws LockedException {
2022 if (running) {
2023 throw new LockedException();
2024 }
2025 this.initialBiasY = convertAngularSpeed(initialBiasY);
2026 }
2027
2028 /**
2029 * Gets initial z-coordinate of gyroscope bias to be used to find a
2030 * solution.
2031 *
2032 * @return initial z-coordinate of gyroscope bias.
2033 */
2034 public AngularSpeed getInitialBiasAngularSpeedZ() {
2035 return new AngularSpeed(initialBiasZ, AngularSpeedUnit.RADIANS_PER_SECOND);
2036 }
2037
2038 /**
2039 * Gets initial z-coordinate of gyroscope bias to be used to find a
2040 * solution.
2041 *
2042 * @param result instance where result data will be stored.
2043 */
2044 public void getInitialBiasAngularSpeedZ(final AngularSpeed result) {
2045 result.setValue(initialBiasZ);
2046 result.setUnit(AngularSpeedUnit.RADIANS_PER_SECOND);
2047 }
2048
2049 /**
2050 * Sets initial z-coordinate of gyroscope bias to be used to find a
2051 * solution.
2052 *
2053 * @param initialBiasZ initial z-coordinate of gyroscope bias.
2054 * @throws LockedException if calibrator is currently running.
2055 */
2056 public void setInitialBiasZ(final AngularSpeed initialBiasZ) throws LockedException {
2057 if (running) {
2058 throw new LockedException();
2059 }
2060 this.initialBiasZ = convertAngularSpeed(initialBiasZ);
2061 }
2062
2063 /**
2064 * Sets initial bias coordinates of gyroscope used to find a solution
2065 * expressed in radians per second (rad/s).
2066 *
2067 * @param initialBiasX initial x-coordinate of gyroscope bias.
2068 * @param initialBiasY initial y-coordinate of gyroscope bias.
2069 * @param initialBiasZ initial z-coordinate of gyroscope bias.
2070 * @throws LockedException if calibrator is currently running.
2071 */
2072 public void setInitialBias(
2073 final double initialBiasX, final double initialBiasY, final double initialBiasZ) throws LockedException {
2074 if (running) {
2075 throw new LockedException();
2076 }
2077 this.initialBiasX = initialBiasX;
2078 this.initialBiasY = initialBiasY;
2079 this.initialBiasZ = initialBiasZ;
2080 }
2081
2082 /**
2083 * Sets initial bias coordinates of gyroscope used to find a solution.
2084 *
2085 * @param initialBiasX initial x-coordinate of gyroscope bias.
2086 * @param initialBiasY initial y-coordinate of gyroscope bias.
2087 * @param initialBiasZ initial z-coordinate of gyroscope bias.
2088 * @throws LockedException if calibrator is currently running.
2089 */
2090 public void setInitialBias(
2091 final AngularSpeed initialBiasX, final AngularSpeed initialBiasY, final AngularSpeed initialBiasZ)
2092 throws LockedException {
2093 if (running) {
2094 throw new LockedException();
2095 }
2096 this.initialBiasX = convertAngularSpeed(initialBiasX);
2097 this.initialBiasY = convertAngularSpeed(initialBiasY);
2098 this.initialBiasZ = convertAngularSpeed(initialBiasZ);
2099 }
2100
2101 /**
2102 * Gets initial x scaling factor of gyroscope.
2103 *
2104 * @return initial x scaling factor of gyroscope.
2105 */
2106 @Override
2107 public double getInitialSx() {
2108 return initialSx;
2109 }
2110
2111 /**
2112 * Sets initial x scaling factor of gyroscope.
2113 *
2114 * @param initialSx initial x scaling factor of gyroscope.
2115 * @throws LockedException if calibrator is currently running.
2116 */
2117 @Override
2118 public void setInitialSx(final double initialSx) throws LockedException {
2119 if (running) {
2120 throw new LockedException();
2121 }
2122 this.initialSx = initialSx;
2123 }
2124
2125 /**
2126 * Gets initial y scaling factor of gyroscope.
2127 *
2128 * @return initial y scaling factor of gyroscope.
2129 */
2130 @Override
2131 public double getInitialSy() {
2132 return initialSy;
2133 }
2134
2135 /**
2136 * Sets initial y scaling factor of gyroscope.
2137 *
2138 * @param initialSy initial y scaling factor of gyroscope.
2139 * @throws LockedException if calibrator is currently running.
2140 */
2141 @Override
2142 public void setInitialSy(final double initialSy) throws LockedException {
2143 if (running) {
2144 throw new LockedException();
2145 }
2146 this.initialSy = initialSy;
2147 }
2148
2149 /**
2150 * Gets initial z scaling factor of gyroscope.
2151 *
2152 * @return initial z scaling factor of gyroscope.
2153 */
2154 @Override
2155 public double getInitialSz() {
2156 return initialSz;
2157 }
2158
2159 /**
2160 * Sets initial z scaling factor of gyroscope.
2161 *
2162 * @param initialSz initial z scaling factor of gyroscope.
2163 * @throws LockedException if calibrator is currently running.
2164 */
2165 @Override
2166 public void setInitialSz(final double initialSz) throws LockedException {
2167 if (running) {
2168 throw new LockedException();
2169 }
2170 this.initialSz = initialSz;
2171 }
2172
2173 /**
2174 * Gets initial x-y cross coupling error of gyroscope.
2175 *
2176 * @return initial x-y cross coupling error of gyroscope.
2177 */
2178 @Override
2179 public double getInitialMxy() {
2180 return initialMxy;
2181 }
2182
2183 /**
2184 * Sets initial x-y cross coupling error of gyroscope.
2185 *
2186 * @param initialMxy initial x-y cross coupling error of gyroscope.
2187 * @throws LockedException if calibrator is currently running.
2188 */
2189 @Override
2190 public void setInitialMxy(final double initialMxy) throws LockedException {
2191 if (running) {
2192 throw new LockedException();
2193 }
2194 this.initialMxy = initialMxy;
2195 }
2196
2197 /**
2198 * Gets initial x-z cross coupling error of gyroscope.
2199 *
2200 * @return initial x-z cross coupling error of gyroscope.
2201 */
2202 @Override
2203 public double getInitialMxz() {
2204 return initialMxz;
2205 }
2206
2207 /**
2208 * Sets initial x-z cross coupling error of gyroscope.
2209 *
2210 * @param initialMxz initial x-z cross coupling error of gyroscope.
2211 * @throws LockedException if calibrator is currently running.
2212 */
2213 @Override
2214 public void setInitialMxz(final double initialMxz) throws LockedException {
2215 if (running) {
2216 throw new LockedException();
2217 }
2218 this.initialMxz = initialMxz;
2219 }
2220
2221 /**
2222 * Gets initial y-x cross coupling error of gyroscope.
2223 *
2224 * @return initial y-x cross coupling error of gyroscope.
2225 */
2226 @Override
2227 public double getInitialMyx() {
2228 return initialMyx;
2229 }
2230
2231 /**
2232 * Sets initial y-x cross coupling error of gyroscope.
2233 *
2234 * @param initialMyx initial y-x cross coupling error of gyroscope.
2235 * @throws LockedException if calibrator is currently running.
2236 */
2237 @Override
2238 public void setInitialMyx(final double initialMyx) throws LockedException {
2239 if (running) {
2240 throw new LockedException();
2241 }
2242 this.initialMyx = initialMyx;
2243 }
2244
2245 /**
2246 * Gets initial y-z cross coupling error of gyroscope.
2247 *
2248 * @return initial y-z cross coupling error of gyroscope.
2249 */
2250 @Override
2251 public double getInitialMyz() {
2252 return initialMyz;
2253 }
2254
2255 /**
2256 * Sets initial y-z cross coupling error of gyroscope.
2257 *
2258 * @param initialMyz initial y-z cross coupling error of gyroscope.
2259 * @throws LockedException if calibrator is currently running.
2260 */
2261 @Override
2262 public void setInitialMyz(final double initialMyz) throws LockedException {
2263 if (running) {
2264 throw new LockedException();
2265 }
2266 this.initialMyz = initialMyz;
2267 }
2268
2269 /**
2270 * Gets initial z-x cross coupling error of gyroscope.
2271 *
2272 * @return initial z-x cross coupling error of gyroscope.
2273 */
2274 @Override
2275 public double getInitialMzx() {
2276 return initialMzx;
2277 }
2278
2279 /**
2280 * Sets initial z-x cross coupling error of gyroscope.
2281 *
2282 * @param initialMzx initial z-x cross coupling error of gyroscope.
2283 * @throws LockedException if calibrator is currently running.
2284 */
2285 @Override
2286 public void setInitialMzx(final double initialMzx) throws LockedException {
2287 if (running) {
2288 throw new LockedException();
2289 }
2290 this.initialMzx = initialMzx;
2291 }
2292
2293 /**
2294 * Gets initial z-y cross coupling error of gyroscope.
2295 *
2296 * @return initial z-y cross coupling error of gyroscope.
2297 */
2298 @Override
2299 public double getInitialMzy() {
2300 return initialMzy;
2301 }
2302
2303 /**
2304 * Sets initial z-y cross coupling error of gyroscope.
2305 *
2306 * @param initialMzy initial z-y cross coupling error of gyroscope.
2307 * @throws LockedException if calibrator is currently running.
2308 */
2309 @Override
2310 public void setInitialMzy(final double initialMzy) throws LockedException {
2311 if (running) {
2312 throw new LockedException();
2313 }
2314 this.initialMzy = initialMzy;
2315 }
2316
2317 /**
2318 * Sets initial scaling factors of gyroscope.
2319 *
2320 * @param initialSx initial x scaling factor of gyroscope.
2321 * @param initialSy initial y scaling factor of gyroscope.
2322 * @param initialSz initial z scaling factor of gyroscope.
2323 * @throws LockedException if calibrator is currently running.
2324 */
2325 @Override
2326 public void setInitialScalingFactors(
2327 final double initialSx, final double initialSy, final double initialSz) throws LockedException {
2328 if (running) {
2329 throw new LockedException();
2330 }
2331 this.initialSx = initialSx;
2332 this.initialSy = initialSy;
2333 this.initialSz = initialSz;
2334 }
2335
2336 /**
2337 * Sets initial cross coupling errors of gyroscope.
2338 *
2339 * @param initialMxy initial x-y cross coupling error of gyroscope.
2340 * @param initialMxz initial x-z cross coupling error of gyroscope.
2341 * @param initialMyx initial y-x cross coupling error of gyroscope.
2342 * @param initialMyz initial y-z cross coupling error of gyroscope.
2343 * @param initialMzx initial z-x cross coupling error of gyroscope.
2344 * @param initialMzy initial z-y cross coupling error of gyroscope.
2345 * @throws LockedException if calibrator is currently running.
2346 */
2347 @Override
2348 public void setInitialCrossCouplingErrors(
2349 final double initialMxy, final double initialMxz, final double initialMyx,
2350 final double initialMyz, final double initialMzx, final double initialMzy) throws LockedException {
2351 if (running) {
2352 throw new LockedException();
2353 }
2354 this.initialMxy = initialMxy;
2355 this.initialMxz = initialMxz;
2356 this.initialMyx = initialMyx;
2357 this.initialMyz = initialMyz;
2358 this.initialMzx = initialMzx;
2359 this.initialMzy = initialMzy;
2360 }
2361
2362 /**
2363 * Sets initial scaling factors and cross coupling errors of
2364 * gyroscope.
2365 *
2366 * @param initialSx initial x scaling factor of gyroscope.
2367 * @param initialSy initial y scaling factor of gyroscope.
2368 * @param initialSz initial z scaling factor of gyroscope.
2369 * @param initialMxy initial x-y cross coupling error of gyroscope.
2370 * @param initialMxz initial x-z cross coupling error of gyroscope.
2371 * @param initialMyx initial y-x cross coupling error of gyroscope.
2372 * @param initialMyz initial y-z cross coupling error of gyroscope.
2373 * @param initialMzx initial z-x cross coupling error of gyroscope.
2374 * @param initialMzy initial z-y cross coupling error of gyroscope.
2375 * @throws LockedException if calibrator is currently running.
2376 */
2377 @Override
2378 public void setInitialScalingFactorsAndCrossCouplingErrors(
2379 final double initialSx, final double initialSy, final double initialSz,
2380 final double initialMxy, final double initialMxz, final double initialMyx,
2381 final double initialMyz, final double initialMzx, final double initialMzy) throws LockedException {
2382 if (running) {
2383 throw new LockedException();
2384 }
2385 setInitialScalingFactors(initialSx, initialSy, initialSz);
2386 setInitialCrossCouplingErrors(initialMxy, initialMxz, initialMyx, initialMyz, initialMzx, initialMzy);
2387 }
2388
2389 /**
2390 * Gets initial gyroscope bias to be used to find a solution as
2391 * an array.
2392 * Array values are expressed in radians per second (rad/s).
2393 *
2394 * @return array containing coordinates of initial gyroscope bias.
2395 */
2396 public double[] getInitialBias() {
2397 final var result = new double[BodyKinematics.COMPONENTS];
2398 getInitialBias(result);
2399 return result;
2400 }
2401
2402 /**
2403 * Gets initial gyroscope bias to be used to find a solution as
2404 * an array.
2405 * Array values are expressed in radians per second (rad/s).
2406 *
2407 * @param result instance where result data will be copied to.
2408 * @throws IllegalArgumentException if provided array does not have length 3.
2409 */
2410 public void getInitialBias(final double[] result) {
2411 if (result.length != BodyKinematics.COMPONENTS) {
2412 throw new IllegalArgumentException();
2413 }
2414 result[0] = initialBiasX;
2415 result[1] = initialBiasY;
2416 result[2] = initialBiasZ;
2417 }
2418
2419 /**
2420 * Sets initial gyroscope bias to be used to find a solution as
2421 * an array.
2422 * Array values are expressed in radians per second (rad/s).
2423 *
2424 * @param initialBias initial bias to find a solution.
2425 * @throws LockedException if calibrator is currently running.
2426 * @throws IllegalArgumentException if provided array does not have length 3.
2427 */
2428 public void setInitialBias(final double[] initialBias) throws LockedException {
2429 if (running) {
2430 throw new LockedException();
2431 }
2432
2433 if (initialBias.length != BodyKinematics.COMPONENTS) {
2434 throw new IllegalArgumentException();
2435 }
2436 initialBiasX = initialBias[0];
2437 initialBiasY = initialBias[1];
2438 initialBiasZ = initialBias[2];
2439 }
2440
2441 /**
2442 * Gets initial gyroscope bias to be used to find a solution as a
2443 * column matrix.
2444 * Values are expressed in radians per second (rad/s).
2445 *
2446 * @return initial gyroscope bias to be used to find a solution as a
2447 * column matrix.
2448 */
2449 public Matrix getInitialBiasAsMatrix() {
2450 Matrix result;
2451 try {
2452 result = new Matrix(BodyKinematics.COMPONENTS, 1);
2453 getInitialBiasAsMatrix(result);
2454 } catch (final WrongSizeException ignore) {
2455 // never happens
2456 result = null;
2457 }
2458 return result;
2459 }
2460
2461 /**
2462 * Gets initial gyroscope bias to be used to find a solution as a
2463 * column matrix.
2464 * Values are expressed in radians per second (rad/s).
2465 *
2466 * @param result instance where result data will be copied to.
2467 * @throws IllegalArgumentException if provided matrix is not 3x1.
2468 */
2469 public void getInitialBiasAsMatrix(final Matrix result) {
2470 if (result.getRows() != BodyKinematics.COMPONENTS || result.getColumns() != 1) {
2471 throw new IllegalArgumentException();
2472 }
2473 result.setElementAtIndex(0, initialBiasX);
2474 result.setElementAtIndex(1, initialBiasY);
2475 result.setElementAtIndex(2, initialBiasZ);
2476 }
2477
2478 /**
2479 * Sets initial gyroscope bias to be used to find a solution as
2480 * a column matrix.
2481 * Values are expressed in radians per second (rad/s).
2482 *
2483 * @param initialBias initial gyroscope bias to find a solution.
2484 * @throws LockedException if calibrator is currently running.
2485 * @throws IllegalArgumentException if provided matrix is not 3x1.
2486 */
2487 public void setInitialBias(final Matrix initialBias) throws LockedException {
2488 if (running) {
2489 throw new LockedException();
2490 }
2491 if (initialBias.getRows() != BodyKinematics.COMPONENTS || initialBias.getColumns() != 1) {
2492 throw new IllegalArgumentException();
2493 }
2494
2495 initialBiasX = initialBias.getElementAtIndex(0);
2496 initialBiasY = initialBias.getElementAtIndex(1);
2497 initialBiasZ = initialBias.getElementAtIndex(2);
2498 }
2499
2500 /**
2501 * Gets initial bias coordinates of gyroscope used to find a solution.
2502 *
2503 * @return initial bias coordinates.
2504 */
2505 public AngularSpeedTriad getInitialBiasAsTriad() {
2506 return new AngularSpeedTriad(AngularSpeedUnit.RADIANS_PER_SECOND, initialBiasX, initialBiasY, initialBiasZ);
2507 }
2508
2509 /**
2510 * Gets initial bias coordinates of gyroscope used to find a solution.
2511 *
2512 * @param result instance where result will be stored.
2513 */
2514 public void getInitialBiasAsTriad(AngularSpeedTriad result) {
2515 result.setValueCoordinatesAndUnit(initialBiasX, initialBiasY, initialBiasZ,
2516 AngularSpeedUnit.RADIANS_PER_SECOND);
2517 }
2518
2519 /**
2520 * Sets initial bias coordinates of gyroscope used to find a solution.
2521 *
2522 * @param initialBias initial bias coordinates to be set.
2523 * @throws LockedException if calibrator is currently running.
2524 */
2525 public void setInitialBias(AngularSpeedTriad initialBias) throws LockedException {
2526 if (running) {
2527 throw new LockedException();
2528 }
2529
2530 initialBiasX = convertAngularSpeed(initialBias.getValueX(), initialBias.getUnit());
2531 initialBiasY = convertAngularSpeed(initialBias.getValueY(), initialBias.getUnit());
2532 initialBiasZ = convertAngularSpeed(initialBias.getValueZ(), initialBias.getUnit());
2533 }
2534
2535 /**
2536 * Gets initial gyroscope scale factors and cross coupling errors
2537 * matrix.
2538 *
2539 * @return initial gyroscope scale factors and cross coupling errors
2540 * matrix.
2541 */
2542 @Override
2543 public Matrix getInitialMg() {
2544 Matrix result;
2545 try {
2546 result = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
2547 getInitialMg(result);
2548 } catch (final WrongSizeException ignore) {
2549 // never happens
2550 result = null;
2551 }
2552 return result;
2553 }
2554
2555 /**
2556 * Gets initial gyroscope scale factors and cross coupling errors
2557 * matrix.
2558 *
2559 * @param result instance where data will be stored.
2560 * @throws IllegalArgumentException if provided matrix is not 3x3.
2561 */
2562 @Override
2563 public void getInitialMg(final Matrix result) {
2564 if (result.getRows() != BodyKinematics.COMPONENTS || result.getColumns() != BodyKinematics.COMPONENTS) {
2565 throw new IllegalArgumentException();
2566 }
2567 result.setElementAtIndex(0, initialSx);
2568 result.setElementAtIndex(1, initialMyx);
2569 result.setElementAtIndex(2, initialMzx);
2570
2571 result.setElementAtIndex(3, initialMxy);
2572 result.setElementAtIndex(4, initialSy);
2573 result.setElementAtIndex(5, initialMzy);
2574
2575 result.setElementAtIndex(6, initialMxz);
2576 result.setElementAtIndex(7, initialMyz);
2577 result.setElementAtIndex(8, initialSz);
2578 }
2579
2580 /**
2581 * Sets initial gyroscope scale factors and cross coupling errors matrix.
2582 *
2583 * @param initialMg initial scale factors and cross coupling errors matrix.
2584 * @throws IllegalArgumentException if provided matrix is not 3x3.
2585 * @throws LockedException if calibrator is currently running.
2586 */
2587 @Override
2588 public void setInitialMg(final Matrix initialMg) throws LockedException {
2589 if (running) {
2590 throw new LockedException();
2591 }
2592 if (initialMg.getRows() != BodyKinematics.COMPONENTS || initialMg.getColumns() != BodyKinematics.COMPONENTS) {
2593 throw new IllegalArgumentException();
2594 }
2595
2596 initialSx = initialMg.getElementAtIndex(0);
2597 initialMyx = initialMg.getElementAtIndex(1);
2598 initialMzx = initialMg.getElementAtIndex(2);
2599
2600 initialMxy = initialMg.getElementAtIndex(3);
2601 initialSy = initialMg.getElementAtIndex(4);
2602 initialMzy = initialMg.getElementAtIndex(5);
2603
2604 initialMxz = initialMg.getElementAtIndex(6);
2605 initialMyz = initialMg.getElementAtIndex(7);
2606 initialSz = initialMg.getElementAtIndex(8);
2607 }
2608
2609 /**
2610 * Gets initial G-dependent cross biases introduced on the gyroscope by the
2611 * specific forces sensed by the accelerometer.
2612 *
2613 * @return a 3x3 matrix containing initial g-dependent cross biases.
2614 */
2615 @Override
2616 public Matrix getInitialGg() {
2617 return new Matrix(initialGg);
2618 }
2619
2620 /**
2621 * Gets initial G-dependent cross biases introduced on the gyroscope by the
2622 * specific forces sensed by the accelerometer.
2623 *
2624 * @param result instance where data will be stored.
2625 * @throws IllegalArgumentException if provided matrix is not 3x3.
2626 */
2627 @Override
2628 public void getInitialGg(final Matrix result) {
2629 if (result.getRows() != BodyKinematics.COMPONENTS || result.getColumns() != BodyKinematics.COMPONENTS) {
2630 throw new IllegalArgumentException();
2631 }
2632
2633 result.copyFrom(initialGg);
2634 }
2635
2636 /**
2637 * Sets initial G-dependent cross biases introduced on the gyroscope by the
2638 * specific forces sensed by the accelerometer.
2639 *
2640 * @param initialGg g-dependent cross biases.
2641 * @throws LockedException if calibrator is currently running.
2642 * @throws IllegalArgumentException if provided matrix is not 3x3.
2643 */
2644 @Override
2645 public void setInitialGg(final Matrix initialGg) throws LockedException {
2646 if (running) {
2647 throw new LockedException();
2648 }
2649
2650 if (initialGg.getRows() != BodyKinematics.COMPONENTS || initialGg.getColumns() != BodyKinematics.COMPONENTS) {
2651 throw new IllegalArgumentException();
2652 }
2653
2654 initialGg.copyTo(this.initialGg);
2655 }
2656
2657 /**
2658 * Gets collection of sequences of timestamped body kinematics
2659 * measurements taken at a given position where the device moves freely
2660 * with different orientations.
2661 *
2662 * @return collection of sequences of timestamped body kinematics
2663 * measurements.
2664 */
2665 @Override
2666 public List<BodyKinematicsSequence<StandardDeviationTimedBodyKinematics>> getSequences() {
2667 return sequences;
2668 }
2669
2670 /**
2671 * Sets collection of sequences of timestamped body kinematics
2672 * measurements taken at a given position where the device moves freely
2673 * with different orientations.
2674 *
2675 * @param sequences collection of sequences of timestamped body
2676 * kinematics measurements.
2677 * @throws LockedException if calibrator is currently running.
2678 */
2679 @Override
2680 public void setSequences(
2681 final List<BodyKinematicsSequence<StandardDeviationTimedBodyKinematics>> sequences) throws LockedException {
2682 if (running) {
2683 throw new LockedException();
2684 }
2685 this.sequences = sequences;
2686 }
2687
2688 /**
2689 * Indicates the type of measurement or sequence used by this calibrator.
2690 *
2691 * @return type of measurement or sequence used by this calibrator.
2692 */
2693 @Override
2694 public GyroscopeCalibratorMeasurementOrSequenceType getMeasurementOrSequenceType() {
2695 return GyroscopeCalibratorMeasurementOrSequenceType.BODY_KINEMATICS_SEQUENCE;
2696 }
2697
2698 /**
2699 * Indicates whether this calibrator requires ordered measurements or sequences
2700 * in a list or not.
2701 *
2702 * @return true if measurements or sequences must be ordered, false otherwise.
2703 */
2704 @Override
2705 public boolean isOrderedMeasurementsOrSequencesRequired() {
2706 return true;
2707 }
2708
2709 /**
2710 * Indicates whether this calibrator requires quality scores for each
2711 * measurement/sequence or not.
2712 *
2713 * @return true if quality scores are required, false otherwise.
2714 */
2715 @Override
2716 public boolean isQualityScoresRequired() {
2717 return false;
2718 }
2719
2720 /**
2721 * Indicates whether z-axis is assumed to be common for accelerometer and
2722 * gyroscope.
2723 * When enabled, this eliminates 3 variables from Ma matrix.
2724 *
2725 * @return true if z-axis is assumed to be common for accelerometer and gyroscope,
2726 * false otherwise.
2727 */
2728 @Override
2729 public boolean isCommonAxisUsed() {
2730 return commonAxisUsed;
2731 }
2732
2733 /**
2734 * Specifies whether z-axis is assumed to be common for accelerometer and
2735 * gyroscope.
2736 * When enabled, this eliminates 3 variables from Ma matrix.
2737 *
2738 * @param commonAxisUsed true if z-axis is assumed to be common for accelerometer
2739 * and gyroscope, false otherwise.
2740 * @throws LockedException if calibrator is currently running.
2741 */
2742 @Override
2743 public void setCommonAxisUsed(final boolean commonAxisUsed) throws LockedException {
2744 if (running) {
2745 throw new LockedException();
2746 }
2747
2748 this.commonAxisUsed = commonAxisUsed;
2749 }
2750
2751 /**
2752 * Indicates whether G-dependent cross biases are being estimated
2753 * or not.
2754 * When enabled, this adds 9 variables from Gg matrix.
2755 *
2756 * @return true if G-dependent cross biases will be estimated,
2757 * false otherwise.
2758 */
2759 public boolean isGDependentCrossBiasesEstimated() {
2760 return estimateGDependentCrossBiases;
2761 }
2762
2763 /**
2764 * Specifies whether G-dependent cross biases are being estimated
2765 * or not.
2766 * When enabled, this adds 9 variables from Gg matrix.
2767 *
2768 * @param estimateGDependentCrossBiases true if G-dependent cross
2769 * biases will be estimated,
2770 * false otherwise.
2771 * @throws LockedException if calibrator is currently running.
2772 */
2773 public void setGDependentCrossBiasesEstimated(final boolean estimateGDependentCrossBiases) throws LockedException {
2774 if (running) {
2775 throw new LockedException();
2776 }
2777
2778 this.estimateGDependentCrossBiases = estimateGDependentCrossBiases;
2779 }
2780
2781 /**
2782 * Gets listener to handle events raised by this estimator.
2783 *
2784 * @return listener to handle events raised by this estimator.
2785 */
2786 public EasyGyroscopeCalibratorListener getListener() {
2787 return listener;
2788 }
2789
2790 /**
2791 * Sets listener to handle events raised by this estimator.
2792 *
2793 * @param listener listener to handle events raised by this estimator.
2794 * @throws LockedException if calibrator is currently running.
2795 */
2796 public void setListener(final EasyGyroscopeCalibratorListener listener) throws LockedException {
2797 if (running) {
2798 throw new LockedException();
2799 }
2800
2801 this.listener = listener;
2802 }
2803
2804 /**
2805 * Gets minimum number of required sequences.
2806 *
2807 * @return minimum number of required sequences.
2808 */
2809 @Override
2810 public int getMinimumRequiredMeasurementsOrSequences() {
2811 if (commonAxisUsed) {
2812 if (estimateGDependentCrossBiases) {
2813 return MINIMUM_SEQUENCES_COMMON_Z_AXIS_AND_CROSS_BIASES;
2814 } else {
2815 return MINIMUM_SEQUENCES_COMMON_Z_AXIS;
2816 }
2817 } else {
2818 if (estimateGDependentCrossBiases) {
2819 return MINIMUM_SEQUENCES_GENERAL_AND_CROSS_BIASES;
2820 } else {
2821 return MINIMUM_SEQUENCES_GENERAL;
2822 }
2823 }
2824 }
2825
2826 /**
2827 * Indicates whether calibrator is ready to start.
2828 *
2829 * @return true if calibrator is ready, false otherwise.
2830 */
2831 @Override
2832 public boolean isReady() {
2833 return sequences != null && sequences.size() >= getMinimumRequiredMeasurementsOrSequences();
2834 }
2835
2836 /**
2837 * Indicates whether calibrator is currently running or not.
2838 *
2839 * @return true if calibrator is running, false otherwise.
2840 */
2841 @Override
2842 public boolean isRunning() {
2843 return running;
2844 }
2845
2846 /**
2847 * Estimates gyroscope calibration parameters containing bias, scale factors,
2848 * cross-coupling errors and G-dependent coupling.
2849 *
2850 * @throws LockedException if calibrator is currently running.
2851 * @throws NotReadyException if calibrator is not ready.
2852 * @throws CalibrationException if estimation fails for numerical reasons.
2853 */
2854 @Override
2855 public void calibrate() throws LockedException, NotReadyException, CalibrationException {
2856 if (running) {
2857 throw new LockedException();
2858 }
2859
2860 if (!isReady()) {
2861 throw new NotReadyException();
2862 }
2863
2864 try {
2865 running = true;
2866
2867 reset();
2868
2869 if (listener != null) {
2870 listener.onCalibrateStart(this);
2871 }
2872
2873 if (commonAxisUsed) {
2874 if (estimateGDependentCrossBiases) {
2875 calibrateCommonAxisAndGDependentCrossBiases();
2876 } else {
2877 calibrateCommonAxis();
2878 }
2879 } else {
2880 if (estimateGDependentCrossBiases) {
2881 calibrateGeneralAndGDependentCrossBiases();
2882 } else {
2883 calibrateGeneral();
2884 }
2885 }
2886
2887 if (listener != null) {
2888 listener.onCalibrateEnd(this);
2889 }
2890
2891 } catch (final FittingException | AlgebraException | com.irurueta.numerical.NotReadyException e) {
2892 throw new CalibrationException(e);
2893 } finally {
2894 running = false;
2895 }
2896 }
2897
2898 /**
2899 * Gets array containing x,y,z components of estimated gyroscope biases
2900 * expressed in radians per second (rad/s).
2901 *
2902 * @return array containing x,y,z components of estimated gyroscope biases.
2903 */
2904 @Override
2905 public double[] getEstimatedBiases() {
2906 return estimatedBiases;
2907 }
2908
2909 /**
2910 * Gets array containing x,y,z components of estimated gyroscope biases
2911 * expressed in radians per second (rad/s).
2912 *
2913 * @param result instance where estimated gyroscope biases will be stored.
2914 * @return true if result instance was updated, false otherwise (when estimation
2915 * is not yet available).
2916 */
2917 @Override
2918 public boolean getEstimatedBiases(final double[] result) {
2919 if (estimatedBiases != null) {
2920 System.arraycopy(estimatedBiases, 0, result, 0, estimatedBiases.length);
2921 return true;
2922 } else {
2923 return false;
2924 }
2925 }
2926
2927 /**
2928 * Gets column matrix containing x,y,z components of estimated gyroscope biases
2929 * expressed in radians per second (rad/s).
2930 *
2931 * @return column matrix containing x,y,z components of estimated gyroscope
2932 * biases.
2933 */
2934 @Override
2935 public Matrix getEstimatedBiasesAsMatrix() {
2936 return estimatedBiases != null ? Matrix.newFromArray(estimatedBiases) : null;
2937 }
2938
2939 /**
2940 * Gets column matrix containing x,y,z components of estimated gyroscope biases
2941 * expressed in radians per second (rad/s).
2942 *
2943 * @param result instance where result data will be stored.
2944 * @return true if result was updated, false otherwise.
2945 * @throws WrongSizeException if provided result instance has invalid size.
2946 */
2947 @Override
2948 public boolean getEstimatedBiasesAsMatrix(final Matrix result) throws WrongSizeException {
2949 if (estimatedBiases != null) {
2950 result.fromArray(estimatedBiases);
2951 return true;
2952 } else {
2953 return false;
2954 }
2955 }
2956
2957 /**
2958 * Gets x coordinate of estimated gyroscope bias expressed in radians per
2959 * second (rad/s).
2960 *
2961 * @return x coordinate of estimated gyroscope bias or null if not available.
2962 */
2963 @Override
2964 public Double getEstimatedBiasX() {
2965 return estimatedBiases != null ? estimatedBiases[0] : null;
2966 }
2967
2968 /**
2969 * Gets y coordinate of estimated gyroscope bias expressed in radians per
2970 * second (rad/s).
2971 *
2972 * @return y coordinate of estimated gyroscope bias or null if not available.
2973 */
2974 @Override
2975 public Double getEstimatedBiasY() {
2976 return estimatedBiases != null ? estimatedBiases[1] : null;
2977 }
2978
2979 /**
2980 * Gets z coordinate of estimated gyroscope bias expressed in radians per
2981 * second (rad/s).
2982 *
2983 * @return z coordinate of estimated gyroscope bias or null if not available.
2984 */
2985 @Override
2986 public Double getEstimatedBiasZ() {
2987 return estimatedBiases != null ? estimatedBiases[2] : null;
2988 }
2989
2990 /**
2991 * Gets x coordinate of estimated gyroscope bias.
2992 *
2993 * @return x coordinate of estimated gyroscope bias or null if not available.
2994 */
2995 @Override
2996 public AngularSpeed getEstimatedBiasAngularSpeedX() {
2997 return estimatedBiases != null
2998 ? new AngularSpeed(estimatedBiases[0], AngularSpeedUnit.RADIANS_PER_SECOND) : null;
2999 }
3000
3001 /**
3002 * Gets x coordinate of estimated gyroscope bias.
3003 *
3004 * @param result instance where result will be stored.
3005 * @return true if result was updated, false if estimation is not available.
3006 */
3007 @Override
3008 public boolean getEstimatedBiasAngularSpeedX(final AngularSpeed result) {
3009 if (estimatedBiases != null) {
3010 result.setValue(estimatedBiases[0]);
3011 result.setUnit(AngularSpeedUnit.RADIANS_PER_SECOND);
3012 return true;
3013 } else {
3014 return false;
3015 }
3016 }
3017
3018 /**
3019 * Gets y coordinate of estimated gyroscope bias.
3020 *
3021 * @return y coordinate of estimated gyroscope bias or null if not available.
3022 */
3023 @Override
3024 public AngularSpeed getEstimatedBiasAngularSpeedY() {
3025 return estimatedBiases != null
3026 ? new AngularSpeed(estimatedBiases[1], AngularSpeedUnit.RADIANS_PER_SECOND) : null;
3027 }
3028
3029 /**
3030 * Gets y coordinate of estimated gyroscope bias.
3031 *
3032 * @param result instance where result will be stored.
3033 * @return true if result was updated, false if estimation is not available.
3034 */
3035 @Override
3036 public boolean getEstimatedBiasAngularSpeedY(final AngularSpeed result) {
3037 if (estimatedBiases != null) {
3038 result.setValue(estimatedBiases[1]);
3039 result.setUnit(AngularSpeedUnit.RADIANS_PER_SECOND);
3040 return true;
3041 } else {
3042 return false;
3043 }
3044 }
3045
3046 /**
3047 * Gets z coordinate of estimated gyroscope bias.
3048 *
3049 * @return z coordinate of estimated gyroscope bias or null if not available.
3050 */
3051 @Override
3052 public AngularSpeed getEstimatedBiasAngularSpeedZ() {
3053 return estimatedBiases != null
3054 ? new AngularSpeed(estimatedBiases[2], AngularSpeedUnit.RADIANS_PER_SECOND) : null;
3055 }
3056
3057 /**
3058 * Gets z coordinate of estimated gyroscope bias.
3059 *
3060 * @param result instance where result will be stored.
3061 * @return true if result was updated, false if estimation is not available.
3062 */
3063 @Override
3064 public boolean getEstimatedBiasAngularSpeedZ(final AngularSpeed result) {
3065 if (estimatedBiases != null) {
3066 result.setValue(estimatedBiases[2]);
3067 result.setUnit(AngularSpeedUnit.RADIANS_PER_SECOND);
3068 return true;
3069 } else {
3070 return false;
3071 }
3072 }
3073
3074 /**
3075 * Gets estimated gyroscope bias.
3076 *
3077 * @return estimated gyroscope bias or null if not available.
3078 */
3079 @Override
3080 public AngularSpeedTriad getEstimatedBiasAsTriad() {
3081 return estimatedBiases != null
3082 ? new AngularSpeedTriad(AngularSpeedUnit.RADIANS_PER_SECOND,
3083 estimatedBiases[0], estimatedBiases[1], estimatedBiases[2])
3084 : null;
3085 }
3086
3087 /**
3088 * Gets estimated gyroscope bias.
3089 *
3090 * @param result instance where result will be stored.
3091 * @return true if estimated gyroscope bias is available and result was
3092 * modified, false otherwise.
3093 */
3094 @Override
3095 public boolean getEstimatedBiasAsTriad(final AngularSpeedTriad result) {
3096 if (estimatedBiases != null) {
3097 result.setValueCoordinatesAndUnit(
3098 estimatedBiases[0], estimatedBiases[1], estimatedBiases[2], AngularSpeedUnit.RADIANS_PER_SECOND);
3099 return true;
3100 } else {
3101 return false;
3102 }
3103 }
3104
3105 /**
3106 * Gets estimated gyroscope scale factors and cross coupling errors.
3107 * This is the product of matrix Tg containing cross coupling errors and Kg
3108 * containing scaling factors.
3109 * So that:
3110 * <pre>
3111 * Mg = [sx mxy mxz] = Tg*Kg
3112 * [myx sy myz]
3113 * [mzx mzy sz ]
3114 * </pre>
3115 * Where:
3116 * <pre>
3117 * Kg = [sx 0 0 ]
3118 * [0 sy 0 ]
3119 * [0 0 sz]
3120 * </pre>
3121 * and
3122 * <pre>
3123 * Tg = [1 -alphaXy alphaXz ]
3124 * [alphaYx 1 -alphaYz]
3125 * [-alphaZx alphaZy 1 ]
3126 * </pre>
3127 * Hence:
3128 * <pre>
3129 * Mg = [sx mxy mxz] = Tg*Kg = [sx -sy * alphaXy sz * alphaXz ]
3130 * [myx sy myz] [sx * alphaYx sy -sz * alphaYz]
3131 * [mzx mzy sz ] [-sx * alphaZx sy * alphaZy sz ]
3132 * </pre>
3133 * This instance allows any 3x3 matrix however, typically alphaYx, alphaZx and alphaZy
3134 * are considered to be zero if the gyroscope z-axis is assumed to be the same
3135 * as the body z-axis. When this is assumed, myx = mzx = mzy = 0 and the Mg matrix
3136 * becomes upper diagonal:
3137 * <pre>
3138 * Mg = [sx mxy mxz]
3139 * [0 sy myz]
3140 * [0 0 sz ]
3141 * </pre>
3142 * Values of this matrix are unit-less.
3143 *
3144 * @return estimated gyroscope scale factors and cross coupling errors, or null
3145 * if not available.
3146 */
3147 @Override
3148 public Matrix getEstimatedMg() {
3149 return estimatedMg;
3150 }
3151
3152 /**
3153 * Gets estimated gyroscope x-axis scale factor.
3154 *
3155 * @return estimated gyroscope x-axis scale factor or null
3156 * if not available.
3157 */
3158 @Override
3159 public Double getEstimatedSx() {
3160 return estimatedMg != null ? estimatedMg.getElementAt(0, 0) : null;
3161 }
3162
3163 /**
3164 * Gets estimated gyroscope y-axis scale factor.
3165 *
3166 * @return estimated gyroscope y-axis scale factor or null
3167 * if not available.
3168 */
3169 @Override
3170 public Double getEstimatedSy() {
3171 return estimatedMg != null ? estimatedMg.getElementAt(1, 1) : null;
3172 }
3173
3174 /**
3175 * Gets estimated gyroscope z-axis scale factor.
3176 *
3177 * @return estimated gyroscope z-axis scale factor or null
3178 * if not available.
3179 */
3180 @Override
3181 public Double getEstimatedSz() {
3182 return estimatedMg != null ? estimatedMg.getElementAt(2, 2) : null;
3183 }
3184
3185 /**
3186 * Gets estimated gyroscope x-y cross-coupling error.
3187 *
3188 * @return estimated gyroscope x-y cross-coupling error or null
3189 * if not available.
3190 */
3191 @Override
3192 public Double getEstimatedMxy() {
3193 return estimatedMg != null ? estimatedMg.getElementAt(0, 1) : null;
3194 }
3195
3196 /**
3197 * Gets estimated gyroscope x-z cross-coupling error.
3198 *
3199 * @return estimated gyroscope x-z cross-coupling error or null
3200 * if not available.
3201 */
3202 @Override
3203 public Double getEstimatedMxz() {
3204 return estimatedMg != null ? estimatedMg.getElementAt(0, 2) : null;
3205 }
3206
3207 /**
3208 * Gets estimated gyroscope y-x cross-coupling error.
3209 *
3210 * @return estimated gyroscope y-x cross-coupling error or null
3211 * if not available.
3212 */
3213 @Override
3214 public Double getEstimatedMyx() {
3215 return estimatedMg != null ? estimatedMg.getElementAt(1, 0) : null;
3216 }
3217
3218 /**
3219 * Gets estimated gyroscope y-z cross-coupling error.
3220 *
3221 * @return estimated gyroscope y-z cross-coupling error or null
3222 * if not available.
3223 */
3224 @Override
3225 public Double getEstimatedMyz() {
3226 return estimatedMg != null ? estimatedMg.getElementAt(1, 2) : null;
3227 }
3228
3229 /**
3230 * Gets estimated gyroscope z-x cross-coupling error.
3231 *
3232 * @return estimated gyroscope z-x cross-coupling error or null
3233 * if not available.
3234 */
3235 @Override
3236 public Double getEstimatedMzx() {
3237 return estimatedMg != null ? estimatedMg.getElementAt(2, 0) : null;
3238 }
3239
3240 /**
3241 * Gets estimated gyroscope z-y cross-coupling error.
3242 *
3243 * @return estimated gyroscope z-y cross-coupling error or null
3244 * if not available.
3245 */
3246 @Override
3247 public Double getEstimatedMzy() {
3248 return estimatedMg != null ? estimatedMg.getElementAt(2, 1) : null;
3249 }
3250
3251 /**
3252 * Gets estimated G-dependent cross biases introduced on the gyroscope by the
3253 * specific forces sensed by the accelerometer.
3254 * This instance allows any 3x3 matrix.
3255 *
3256 * @return estimated G-dependent cross biases.
3257 */
3258 @Override
3259 public Matrix getEstimatedGg() {
3260 return estimatedGg;
3261 }
3262
3263 /**
3264 * Gets estimated covariance matrix for estimated parameters.
3265 * Diagonal elements of the matrix contains variance for the following
3266 * parameters (following indicated order): bgx, bgy, bgz, sx, sy, sz,
3267 * mxy, mxz, myx, myz, mzx, mzy, gg11, gg21, gg31, gg12, gg22, gg32,
3268 * gg13, gg23, gg33.
3269 *
3270 * @return estimated covariance matrix for estimated parameters.
3271 */
3272 @Override
3273 public Matrix getEstimatedCovariance() {
3274 return estimatedCovariance;
3275 }
3276
3277 /**
3278 * Gets estimated chi square value.
3279 *
3280 * @return estimated chi square value.
3281 */
3282 @Override
3283 public double getEstimatedChiSq() {
3284 return estimatedChiSq;
3285 }
3286
3287 /**
3288 * Gets estimated chi square degrees of freedom. Degrees of freedom is equal to the number of sampled data minus the
3289 * number of estimated parameters.
3290 *
3291 * @return estimated degrees of freedom of chi square value
3292 */
3293 @Override
3294 public int getEstimatedChiSqDegreesOfFreedom() {
3295 return estimatedChiSqDegreesOfFreedom;
3296 }
3297
3298 /**
3299 * Gets estimated reduced chi square value. This is equal to estimated chi square value divided by its degrees of
3300 * freedom. Ideally this value should be close to 1.0, indicating that fit is optimal.
3301 * A value larger than 1.0 indicates that fit is not good or noise has been underestimated, and a value smaller than
3302 * 1.0 indicates that there is overfitting or noise has been overestimated.
3303 *
3304 * @return estimated reduced chi square value
3305 */
3306 @Override
3307 public double getEstimatedReducedChiSq() {
3308 return estimatedReducedChiSq;
3309 }
3310
3311 /**
3312 * Gets estimated probability of finding a smaller chi square value expressed as a value between 0.0 and 1.0. The
3313 * smaller the found chi square value is, the better the fit of the estimated parameters to the actual parameter.
3314 * Thus, the smaller the chance of finding a smaller chi square value, then the better the estimated fit is.
3315 *
3316 * @return estimated probability of finding a smaller chi square value.
3317 */
3318 @Override
3319 public double getEstimatedP() {
3320 return estimatedP;
3321 }
3322
3323 /**
3324 * Gets estimated measure of quality of estimated fit as a value between 0.0 and 1.0. The larger the quality value
3325 * is, the better the fit that has been estimated.
3326 *
3327 * @return estimated measure of quality of estimated fit.
3328 */
3329 @Override
3330 public double getEstimatedQ() {
3331 return estimatedQ;
3332 }
3333
3334 /**
3335 * Gets estimated mean square error respect to provided measurements.
3336 *
3337 * @return estimated mean square error respect to provided measurements.
3338 */
3339 @Override
3340 public double getEstimatedMse() {
3341 return estimatedMse;
3342 }
3343
3344 /**
3345 * Gets variance of estimated x coordinate of gyroscope bias expressed in (rad^2/s^2).
3346 *
3347 * @return variance of estimated x coordinate of gyroscope bias or null if not available.
3348 */
3349 public Double getEstimatedBiasXVariance() {
3350 return estimatedCovariance != null ? estimatedCovariance.getElementAt(0, 0) : null;
3351 }
3352
3353 /**
3354 * Gets standard deviation of estimated x coordinate of gyroscope bias expressed in
3355 * radians per second (rad/s).
3356 *
3357 * @return standard deviation of estimated x coordinate of gyroscope bias or null if not
3358 * available.
3359 */
3360 public Double getEstimatedBiasXStandardDeviation() {
3361 final var variance = getEstimatedBiasXVariance();
3362 return variance != null ? Math.sqrt(variance) : null;
3363 }
3364
3365 /**
3366 * Gets standard deviation of estimated x coordinate of gyroscope bias.
3367 *
3368 * @return standard deviation of estimated x coordinate of gyroscope bias or null if not
3369 * available.
3370 */
3371 public AngularSpeed getEstimatedBiasXStandardDeviationAsAngularSpeed() {
3372 return estimatedCovariance != null
3373 ? new AngularSpeed(getEstimatedBiasXStandardDeviation(), AngularSpeedUnit.RADIANS_PER_SECOND)
3374 : null;
3375 }
3376
3377 /**
3378 * Gets standard deviation of estimated x coordinate of gyroscope bias.
3379 *
3380 * @param result instance where result will be stored.
3381 * @return true if standard deviation of estimated x coordinate of gyroscope bias is available,
3382 * false otherwise.
3383 */
3384 public boolean getEstimatedBiasXStandardDeviationAsAngularSpeed(final AngularSpeed result) {
3385 if (estimatedCovariance != null) {
3386 result.setValue(getEstimatedBiasXStandardDeviation());
3387 result.setUnit(AngularSpeedUnit.RADIANS_PER_SECOND);
3388 return true;
3389 } else {
3390 return false;
3391 }
3392 }
3393
3394 /**
3395 * Gets variance of estimated y coordinate of gyroscope bias expressed in (rad^2/s^2).
3396 *
3397 * @return variance of estimated y coordinate of gyroscope bias or null if not available.
3398 */
3399 public Double getEstimatedBiasYVariance() {
3400 return estimatedCovariance != null ? estimatedCovariance.getElementAt(1, 1) : null;
3401 }
3402
3403 /**
3404 * Gets standard deviation of estimated y coordinate of gyroscope bias expressed in
3405 * radians per second (rad/s).
3406 *
3407 * @return standard deviation of estimated y coordinate of gyroscope bias or null if not
3408 * available.
3409 */
3410 public Double getEstimatedBiasYStandardDeviation() {
3411 final var variance = getEstimatedBiasYVariance();
3412 return variance != null ? Math.sqrt(variance) : null;
3413 }
3414
3415 /**
3416 * Gets standard deviation of estimated y coordinate of gyroscope bias.
3417 *
3418 * @return standard deviation of estimated y coordinate of gyroscope bias or null if not
3419 * available.
3420 */
3421 public AngularSpeed getEstimatedBiasYStandardDeviationAsAngularSpeed() {
3422 return estimatedCovariance != null
3423 ? new AngularSpeed(getEstimatedBiasYStandardDeviation(), AngularSpeedUnit.RADIANS_PER_SECOND)
3424 : null;
3425 }
3426
3427 /**
3428 * Gets standard deviation of estimated y coordinate of gyroscope bias.
3429 *
3430 * @param result instance where result will be stored.
3431 * @return true if standard deviation of estimated y coordinate of gyroscope bias is available,
3432 * false otherwise.
3433 */
3434 public boolean getEstimatedBiasYStandardDeviationAsAngularSpeed(final AngularSpeed result) {
3435 if (estimatedCovariance != null) {
3436 result.setValue(getEstimatedBiasYStandardDeviation());
3437 result.setUnit(AngularSpeedUnit.RADIANS_PER_SECOND);
3438 return true;
3439 } else {
3440 return false;
3441 }
3442 }
3443
3444 /**
3445 * Gets variance of estimated z coordinate of gyroscope bias expressed in (rad^2/s^2).
3446 *
3447 * @return variance of estimated z coordinate of gyroscope bias or null if not available.
3448 */
3449 public Double getEstimatedBiasZVariance() {
3450 return estimatedCovariance != null ? estimatedCovariance.getElementAt(2, 2) : null;
3451 }
3452
3453 /**
3454 * Gets standard deviation of estimated z coordinate of gyroscope bias expressed in
3455 * radians per second (rad/s).
3456 *
3457 * @return standard deviation of estimated z coordinate of gyroscope bias or null if not
3458 * available.
3459 */
3460 public Double getEstimatedBiasZStandardDeviation() {
3461 final var variance = getEstimatedBiasZVariance();
3462 return variance != null ? Math.sqrt(variance) : null;
3463 }
3464
3465 /**
3466 * Gets standard deviation of estimated z coordinate of gyroscope bias.
3467 *
3468 * @return standard deviation of estimated z coordinate of gyroscope bias or null if not
3469 * available.
3470 */
3471 public AngularSpeed getEstimatedBiasZStandardDeviationAsAngularSpeed() {
3472 return estimatedCovariance != null
3473 ? new AngularSpeed(getEstimatedBiasZStandardDeviation(), AngularSpeedUnit.RADIANS_PER_SECOND)
3474 : null;
3475 }
3476
3477 /**
3478 * Gets standard deviation of estimated z coordinate of gyroscope bias.
3479 *
3480 * @param result instance where result will be stored.
3481 * @return true if standard deviation of estimated z coordinate of gyroscope bias is available,
3482 * false otherwise.
3483 */
3484 public boolean getEstimatedBiasZStandardDeviationAsAngularSpeed(final AngularSpeed result) {
3485 if (estimatedCovariance != null) {
3486 result.setValue(getEstimatedBiasZStandardDeviation());
3487 result.setUnit(AngularSpeedUnit.RADIANS_PER_SECOND);
3488 return true;
3489 } else {
3490 return false;
3491 }
3492 }
3493
3494 /**
3495 * Gets standard deviation of estimated gyroscope bias coordinates.
3496 *
3497 * @return standard deviation of estimated gyroscope bias coordinates.
3498 */
3499 public AngularSpeedTriad getEstimatedBiasStandardDeviation() {
3500 return estimatedCovariance != null
3501 ? new AngularSpeedTriad(AngularSpeedUnit.RADIANS_PER_SECOND,
3502 getEstimatedBiasXStandardDeviation(),
3503 getEstimatedBiasYStandardDeviation(),
3504 getEstimatedBiasZStandardDeviation())
3505 : null;
3506 }
3507
3508 /**
3509 * Gets standard deviation of estimated gyroscope bias coordinates.
3510 *
3511 * @param result instance where result will be stored.
3512 * @return true if standard deviation of gyroscope bias was available, false
3513 * otherwise.
3514 */
3515 public boolean getEstimatedBiasStandardDeviation(final AngularSpeedTriad result) {
3516 if (estimatedCovariance != null) {
3517 result.setValueCoordinatesAndUnit(
3518 getEstimatedBiasXStandardDeviation(),
3519 getEstimatedBiasYStandardDeviation(),
3520 getEstimatedBiasZStandardDeviation(),
3521 AngularSpeedUnit.RADIANS_PER_SECOND);
3522 return true;
3523 } else {
3524 return false;
3525 }
3526 }
3527
3528 /**
3529 * Gets average of estimated standard deviation of gyroscope bias coordinates expressed
3530 * in radians per second (rad/s).
3531 *
3532 * @return average of estimated standard deviation of gyroscope bias coordinates or null
3533 * if not available.
3534 */
3535 public Double getEstimatedBiasStandardDeviationAverage() {
3536 return estimatedCovariance != null
3537 ? (getEstimatedBiasXStandardDeviation() + getEstimatedBiasYStandardDeviation()
3538 + getEstimatedBiasZStandardDeviation()) / 3.0 : null;
3539 }
3540
3541 /**
3542 * Gets average of estimated standard deviation of gyroscope bias coordinates.
3543 *
3544 * @return average of estimated standard deviation of gyroscope bias coordinates or null.
3545 */
3546 public AngularSpeed getEstimatedBiasStandardDeviationAverageAsAngularSpeed() {
3547 return estimatedCovariance != null
3548 ? new AngularSpeed(getEstimatedBiasStandardDeviationAverage(), AngularSpeedUnit.RADIANS_PER_SECOND)
3549 : null;
3550 }
3551
3552 /**
3553 * Gets average of estimated standard deviation of gyroscope bias coordinates.
3554 *
3555 * @param result instance where result will be stored.
3556 * @return true if average of estimated standard deviation of gyroscope bias is available,
3557 * false otherwise.
3558 */
3559 public boolean getEstimatedBiasStandardDeviationAverageAsAngularSpeed(final AngularSpeed result) {
3560 if (estimatedCovariance != null) {
3561 result.setValue(getEstimatedBiasStandardDeviationAverage());
3562 result.setUnit(AngularSpeedUnit.RADIANS_PER_SECOND);
3563 return true;
3564 } else {
3565 return false;
3566 }
3567 }
3568
3569 /**
3570 * Gets norm of estimated standard deviation of gyroscope bias expressed in
3571 * radians per second (rad/s).
3572 * This can be used as the initial gyroscope bias uncertainty for
3573 * {@link INSLooselyCoupledKalmanInitializerConfig} or {@link INSTightlyCoupledKalmanInitializerConfig}.
3574 *
3575 * @return norm of estimated standard deviation of gyroscope bias or null
3576 * if not available.
3577 */
3578 @Override
3579 public Double getEstimatedBiasStandardDeviationNorm() {
3580 return estimatedCovariance != null
3581 ? Math.sqrt(getEstimatedBiasXVariance() + getEstimatedBiasYVariance() + getEstimatedBiasZVariance())
3582 : null;
3583 }
3584
3585 /**
3586 * Gets norm of estimated standard deviation of gyroscope bias.
3587 * This can be used as the initial gyroscope bias uncertainty for
3588 * {@link INSLooselyCoupledKalmanInitializerConfig} or {@link INSTightlyCoupledKalmanInitializerConfig}.
3589 *
3590 * @return norm of estimated standard deviation of gyroscope bias or null
3591 * if not available.
3592 */
3593 public AngularSpeed getEstimatedBiasStandardDeviationNormAsAngularSpeed() {
3594 return estimatedCovariance != null
3595 ? new AngularSpeed(getEstimatedBiasStandardDeviationNorm(), AngularSpeedUnit.RADIANS_PER_SECOND)
3596 : null;
3597 }
3598
3599 /**
3600 * Gets norm of estimated standard deviation of gyroscope bias coordinates.
3601 * This can be used as the initial gyroscope bias uncertainty for
3602 * {@link INSLooselyCoupledKalmanInitializerConfig} or {@link INSTightlyCoupledKalmanInitializerConfig}.
3603 *
3604 * @param result instance where result will be stored.
3605 * @return true if norm of estimated standard deviation of gyroscope bias is
3606 * available, false otherwise.
3607 */
3608 public boolean getEstimatedBiasStandardDeviationNormAsAngularSpeed(final AngularSpeed result) {
3609 if (estimatedCovariance != null) {
3610 result.setValue(getEstimatedBiasStandardDeviationNorm());
3611 result.setUnit(AngularSpeedUnit.RADIANS_PER_SECOND);
3612 return true;
3613 } else {
3614 return false;
3615 }
3616 }
3617
3618 /**
3619 * Internal method to perform calibration when common z-axis is assumed
3620 * for both the accelerometer and gyroscope and when G-dependent cross
3621 * biases are being estimated.
3622 *
3623 * @throws AlgebraException if accelerometer parameters prevent fixing
3624 * measured accelerometer values due to numerical instabilities.
3625 * @throws FittingException if no convergence to solution is found.
3626 * @throws com.irurueta.numerical.NotReadyException if fitter is not ready.
3627 */
3628 private void calibrateCommonAxisAndGDependentCrossBiases() throws AlgebraException, FittingException,
3629 com.irurueta.numerical.NotReadyException {
3630 // The gyroscope model is
3631 // Ωmeas = bg + (I + Mg) * Ωtrue + Gg * ftrue + w
3632
3633 // Ideally a least squares solution tries to minimize noise component, so:
3634 // Ωmeas = bg + (I + Mg) * Ωtrue + Gg * ftrue
3635
3636 // Where Mg is upper triangular
3637
3638 // For convergence purposes of the Levenberg-Marquardt algorithm, we
3639 // take common factor M = I + Mg
3640
3641 // and the gyroscope model can be better expressed as:
3642
3643 // Ωmeas = M*(Ωtrue + b + G * ftrue)
3644
3645 // where:
3646 // bg = M*b --> b = M^-1*bg
3647 // Gg = M*G --> G = M^-1*Gg
3648
3649 // Hence
3650 // Ωmeas - M*b - M*G*ftrue = M*Ωtrue
3651 // M^-1 * (Ωmeas - M*b - M*G*ftrue) = Ωtrue
3652
3653 // Notice that M is upper diagonal because Mg is upper diagonal
3654 // when common axis is assumed
3655
3656 final var gradientEstimator = new GradientEstimator(
3657 params -> evaluateCommonAxisWithGDependentCrossBiases(i, params));
3658
3659 final var initialM = Matrix.identity(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
3660 initialM.add(getInitialMg());
3661
3662 // Force initial M to be upper diagonal
3663 initialM.setElementAt(1, 0, 0.0);
3664 initialM.setElementAt(2, 0, 0.0);
3665 initialM.setElementAt(2, 1, 0.0);
3666
3667 final var invInitialM = Utils.inverse(initialM);
3668 final var initialBg = getInitialBiasAsMatrix();
3669 final var initialB = invInitialM.multiplyAndReturnNew(initialBg);
3670 final var initialGg = getInitialGg();
3671 final var initialG = invInitialM.multiplyAndReturnNew(initialGg);
3672
3673 fitter.setFunctionEvaluator(new LevenbergMarquardtMultiDimensionFunctionEvaluator() {
3674 @Override
3675 public int getNumberOfDimensions() {
3676 // Before and after normalized gravity versors
3677 return 2 * BodyKinematics.COMPONENTS;
3678 }
3679
3680 @Override
3681 public double[] createInitialParametersArray() {
3682 final var initial = new double[COMMON_Z_AXIS_UNKNOWNS_AND_CROSS_BIASES];
3683
3684 // biases b
3685 for (var i = 0; i < BodyKinematics.COMPONENTS; i++) {
3686 initial[i] = initialB.getElementAtIndex(i);
3687 }
3688
3689 // upper diagonal cross coupling errors M
3690 var k = BodyKinematics.COMPONENTS;
3691 for (var j = 0; j < BodyKinematics.COMPONENTS; j++) {
3692 for (var i = 0; i < BodyKinematics.COMPONENTS; i++) {
3693 if (i <= j) {
3694 initial[k] = initialM.getElementAt(i, j);
3695 k++;
3696 }
3697 }
3698 }
3699
3700 // g-dependent cross biases G
3701 final var num = BodyKinematics.COMPONENTS * BodyKinematics.COMPONENTS;
3702 for (int i = 0, j = k; i < num; i++, j++) {
3703 initial[j] = initialG.getElementAtIndex(i);
3704 }
3705
3706 return initial;
3707 }
3708
3709 @Override
3710 public double evaluate(
3711 final int i, final double[] point, final double[] params, final double[] derivatives)
3712 throws EvaluationException {
3713 EasyGyroscopeCalibrator.this.i = i;
3714
3715 // point contains fixed gravity versor values for current
3716 // sequence
3717 EasyGyroscopeCalibrator.this.point = point;
3718
3719 gradientEstimator.gradient(params, derivatives);
3720
3721 return evaluateCommonAxisWithGDependentCrossBiases(i, params);
3722 }
3723 });
3724
3725 setInputData();
3726
3727 fitter.fit();
3728
3729 final var result = fitter.getA();
3730
3731 final var bx = result[0];
3732 final var by = result[1];
3733 final var bz = result[2];
3734
3735 final var m11 = result[3];
3736
3737 final var m12 = result[4];
3738 final var m22 = result[5];
3739
3740 final var m13 = result[6];
3741 final var m23 = result[7];
3742 final var m33 = result[8];
3743
3744 final var g11 = result[9];
3745 final var g21 = result[10];
3746 final var g31 = result[11];
3747
3748 final var g12 = result[12];
3749 final var g22 = result[13];
3750 final var g32 = result[14];
3751
3752 final var g13 = result[15];
3753 final var g23 = result[16];
3754 final var g33 = result[17];
3755
3756 final var b = new Matrix(BodyKinematics.COMPONENTS, 1);
3757 b.setElementAtIndex(0, bx);
3758 b.setElementAtIndex(1, by);
3759 b.setElementAtIndex(2, bz);
3760
3761 final var m = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
3762 m.setElementAtIndex(0, m11);
3763 m.setElementAtIndex(1, 0.0);
3764 m.setElementAtIndex(2, 0.0);
3765
3766 m.setElementAtIndex(3, m12);
3767 m.setElementAtIndex(4, m22);
3768 m.setElementAtIndex(5, 0.0);
3769
3770 m.setElementAtIndex(6, m13);
3771 m.setElementAtIndex(7, m23);
3772 m.setElementAtIndex(8, m33);
3773
3774 final var g = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
3775 g.setElementAtIndex(0, g11);
3776 g.setElementAtIndex(1, g21);
3777 g.setElementAtIndex(2, g31);
3778
3779 g.setElementAtIndex(3, g12);
3780 g.setElementAtIndex(4, g22);
3781 g.setElementAtIndex(5, g32);
3782
3783 g.setElementAtIndex(6, g13);
3784 g.setElementAtIndex(7, g23);
3785 g.setElementAtIndex(8, g33);
3786
3787 setResult(m, b, g);
3788
3789 // at this point covariance is expressed in terms of b, M and G, and must
3790 // be expressed in terms of bg, Mg and Gg.
3791 // We know that:
3792 // bg = M * b
3793 // Mg = M - I
3794 // Gg = M * G
3795
3796 // b = [bx]
3797 // [by]
3798 // [bz]
3799
3800 // M = [m11 m12 m13]
3801 // [0 m22 m23]
3802 // [0 0 m33]
3803
3804 // G = [g11 g12 g13]
3805 // [g21 g22 g23]
3806 // [g31 g32 g33]
3807
3808 // bg = [m11 m12 m13][bx] = [m11 * bx + m12 * by + m13 * bz] = [bgx]
3809 // [0 m22 m23][by] [ m22 * by + m23 * bz] [bgy]
3810 // [0 0 m33][bz] [ m33 * bz] [bgz]
3811
3812 // Mg = [sx mxy mxz] = [m11 - 1 m12 m13 ]
3813 // [myx sy myz] [0 m22 - 1 m23 ]
3814 // [mzx mzy sz ] [0 0 m33 - 1 ]
3815
3816 // Gg = [gg11 gg12 gg13] = [m11 m12 m13][g11 g12 g13]
3817 // [gg21 gg22 gg23] [0 m22 m23][g21 g22 g23]
3818 // [gg31 gg32 gg33] [0 0 m33][g31 g32 g33]
3819
3820 // Defining the linear application:
3821 // F(b, M, G) = F(bx, by, bz, m11, m12, m22, m13, m23, m33, g11, g21, g31, g12, g22, g32, g13, g23, g33)
3822 // as:
3823 // [bgx] = [m11 * bx + m12 * by + m13 * bz]
3824 // [bgy] [ m22 * by + m23 * bz]
3825 // [bgz] [ m33 * bz]
3826 // [sx] [m11 - 1]
3827 // [sy] [m22 - 1]
3828 // [sz] [m33 - 1]
3829 // [mxy] [m12]
3830 // [mxz] [m13]
3831 // [myx] [0]
3832 // [myz] [m23]
3833 // [mzx] [0]
3834 // [mzy] [0]
3835 // [gg11] [m11 * g11 + m12 * g21 + m13 * g31]
3836 // [gg21] [ m22 * g21 + m23 * g31]
3837 // [gg31] [ m33 * g31]
3838 // [gg12] [m11 * g12 + m12 * g22 + m13 * g32]
3839 // [gg22] [ m22 * g22 + m23 * g32]
3840 // [gg32] [ m33 * g32]
3841 // [gg13] [m11 * g13 + m12 * g23 + m13 * g33]
3842 // [gg23] [ m22 * g23 + m23 * g33]
3843 // [gg33] [ m33 * g33]
3844
3845 // Then the Jacobian of F(b, M, G) is:
3846 // J = [m11 m12 m13 bx by 0 bz 0 0 0 0 0 0 0 0 0 0 0 ]
3847 // [0 m22 m23 0 0 by 0 bz 0 0 0 0 0 0 0 0 0 0 ]
3848 // [0 0 m33 0 0 0 0 0 bz 0 0 0 0 0 0 0 0 0 ]
3849 // [0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ]
3850 // [0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 ]
3851 // [0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 ]
3852 // [0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 ]
3853 // [0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 ]
3854 // [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ]
3855 // [0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 ]
3856 // [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ]
3857 // [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ]
3858 // [0 0 0 g11 g21 0 g31 0 0 m11 m12 m13 0 0 0 0 0 0 ]
3859 // [0 0 0 0 0 g21 0 g31 0 0 m22 m23 0 0 0 0 0 0 ]
3860 // [0 0 0 0 0 0 0 0 g31 0 0 m33 0 0 0 0 0 0 ]
3861 // [0 0 0 g12 g22 0 g32 0 0 0 0 0 m11 m12 m13 0 0 0 ]
3862 // [0 0 0 0 0 g22 0 g32 0 0 0 0 0 m22 m23 0 0 0 ]
3863 // [0 0 0 0 0 0 0 0 g32 0 0 0 0 0 m33 0 0 0 ]
3864 // [0 0 0 g13 g23 0 g33 0 0 0 0 0 0 0 0 m11 m12 m13]
3865 // [0 0 0 0 0 g23 0 g33 0 0 0 0 0 0 0 0 m22 m23]
3866 // [0 0 0 0 0 0 0 0 g33 0 0 0 0 0 0 0 0 m33]
3867
3868 // We know that the propagated covariance is J * Cov * J', hence:
3869 final var jacobian = new Matrix(GENERAL_UNKNOWNS_AND_CROSS_BIASES, COMMON_Z_AXIS_UNKNOWNS_AND_CROSS_BIASES);
3870
3871 jacobian.setElementAt(0, 0, m11);
3872 jacobian.setElementAt(0, 1, m12);
3873 jacobian.setElementAt(0, 2, m13);
3874 jacobian.setElementAt(0, 3, bx);
3875 jacobian.setElementAt(0, 4, by);
3876 jacobian.setElementAt(0, 6, bz);
3877
3878 jacobian.setElementAt(1, 1, m22);
3879 jacobian.setElementAt(1, 2, m23);
3880 jacobian.setElementAt(1, 5, by);
3881 jacobian.setElementAt(1, 7, bz);
3882
3883 jacobian.setElementAt(2, 2, m33);
3884 jacobian.setElementAt(2, 8, bz);
3885
3886 jacobian.setElementAt(3, 3, 1.0);
3887 jacobian.setElementAt(4, 5, 1.0);
3888 jacobian.setElementAt(5, 8, 1.0);
3889
3890 jacobian.setElementAt(6, 4, 1.0);
3891 jacobian.setElementAt(7, 6, 1.0);
3892
3893 jacobian.setElementAt(9, 7, 1.0);
3894
3895 jacobian.setElementAt(12, 3, g11);
3896 jacobian.setElementAt(12, 4, g21);
3897 jacobian.setElementAt(12, 6, g31);
3898 jacobian.setElementAt(12, 9, m11);
3899 jacobian.setElementAt(12, 10, m12);
3900 jacobian.setElementAt(12, 11, m13);
3901
3902 jacobian.setElementAt(13, 5, g21);
3903 jacobian.setElementAt(13, 7, g31);
3904 jacobian.setElementAt(13, 10, m22);
3905 jacobian.setElementAt(13, 11, m23);
3906
3907 jacobian.setElementAt(14, 8, g31);
3908 jacobian.setElementAt(14, 11, m33);
3909
3910 jacobian.setElementAt(15, 3, g12);
3911 jacobian.setElementAt(15, 4, g22);
3912 jacobian.setElementAt(15, 6, g32);
3913 jacobian.setElementAt(15, 12, m11);
3914 jacobian.setElementAt(15, 13, m12);
3915 jacobian.setElementAt(15, 14, m13);
3916
3917 jacobian.setElementAt(16, 5, g22);
3918 jacobian.setElementAt(16, 7, g32);
3919 jacobian.setElementAt(16, 13, m22);
3920 jacobian.setElementAt(16, 14, m23);
3921
3922 jacobian.setElementAt(17, 8, g32);
3923 jacobian.setElementAt(17, 14, m33);
3924
3925 jacobian.setElementAt(18, 3, g13);
3926 jacobian.setElementAt(18, 4, g23);
3927 jacobian.setElementAt(18, 6, g33);
3928 jacobian.setElementAt(18, 15, m11);
3929 jacobian.setElementAt(18, 16, m12);
3930 jacobian.setElementAt(18, 17, m13);
3931
3932 jacobian.setElementAt(19, 5, g23);
3933 jacobian.setElementAt(19, 7, g33);
3934 jacobian.setElementAt(19, 16, m22);
3935 jacobian.setElementAt(19, 17, m23);
3936
3937 jacobian.setElementAt(20, 8, g33);
3938 jacobian.setElementAt(20, 17, m33);
3939
3940 final var jacobianTrans = jacobian.transposeAndReturnNew();
3941 jacobian.multiply(estimatedCovariance);
3942 jacobian.multiply(jacobianTrans);
3943 estimatedCovariance = jacobian;
3944 }
3945
3946 /**
3947 * Internal method to perform general calibration when G-dependent cross
3948 * biases are being estimated.
3949 *
3950 * @throws AlgebraException if accelerometer parameters prevent fixing
3951 * measured accelerometer values due to numerical instabilities.
3952 * @throws FittingException if no convergence to solution is found.
3953 * @throws com.irurueta.numerical.NotReadyException if fitter is not ready.
3954 */
3955 private void calibrateGeneralAndGDependentCrossBiases() throws AlgebraException, FittingException,
3956 com.irurueta.numerical.NotReadyException {
3957 // The gyroscope model is
3958 // Ωmeas = bg + (I + Mg) * Ωtrue + Gg * ftrue + w
3959
3960 // Ideally a least squares solution tries to minimize noise component, so:
3961 // Ωmeas = bg + (I + Mg) * Ωtrue + Gg * ftrue
3962
3963 // For convergence purposes of the Levenberg-Marquardt algorithm, we
3964 // take common factor M = I + Mg
3965
3966 // and the gyroscope model can be better expressed as:
3967
3968 // Ωmeas = M*(Ωtrue + b + G * ftrue)
3969
3970 // where:
3971 // bg = M*b --> b = M^-1*bg
3972 // Gg = M*G --> G = M^-1*Gg
3973
3974 // Hence
3975 // Ωmeas - M*b - M*G*ftrue = M*Ωtrue
3976 // M^-1 * (Ωmeas - M*b - M*G*ftrue) = Ωtrue
3977
3978 final var gradientEstimator = new GradientEstimator(
3979 params -> evaluateGeneralWithGDependentCrossBiases(i, params));
3980
3981 final var initialM = Matrix.identity(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
3982 initialM.add(getInitialMg());
3983
3984 final var invInitialM = Utils.inverse(initialM);
3985 final var initialBg = getInitialBiasAsMatrix();
3986 final var initialB = invInitialM.multiplyAndReturnNew(initialBg);
3987 final var initialGg = getInitialGg();
3988 final var initialG = invInitialM.multiplyAndReturnNew(initialGg);
3989
3990 fitter.setFunctionEvaluator(new LevenbergMarquardtMultiDimensionFunctionEvaluator() {
3991 @Override
3992 public int getNumberOfDimensions() {
3993 // Before and after normalized gravity versors
3994 return 2 * BodyKinematics.COMPONENTS;
3995 }
3996
3997 @Override
3998 public double[] createInitialParametersArray() {
3999 final var initial = new double[GENERAL_UNKNOWNS_AND_CROSS_BIASES];
4000
4001 // biases b
4002 for (var i = 0; i < BodyKinematics.COMPONENTS; i++) {
4003 initial[i] = initialB.getElementAtIndex(i);
4004 }
4005
4006 // cross coupling errors M
4007 final var num = BodyKinematics.COMPONENTS * BodyKinematics.COMPONENTS;
4008 for (int i = 0, j = BodyKinematics.COMPONENTS; i < num; i++, j++) {
4009 initial[j] = initialM.getElementAtIndex(i);
4010 }
4011
4012 // g-dependent cross biases G
4013 for (int i = 0, j = BodyKinematics.COMPONENTS + num; i < num; i++, j++) {
4014 initial[j] = initialG.getElementAtIndex(i);
4015 }
4016
4017 return initial;
4018 }
4019
4020 @Override
4021 public double evaluate(
4022 final int i, final double[] point, final double[] params, final double[] derivatives)
4023 throws EvaluationException {
4024 EasyGyroscopeCalibrator.this.i = i;
4025
4026 // point contains fixed gravity versor values for current
4027 // sequence
4028 EasyGyroscopeCalibrator.this.point = point;
4029
4030 gradientEstimator.gradient(params, derivatives);
4031
4032 return evaluateGeneralWithGDependentCrossBiases(i, params);
4033 }
4034 });
4035
4036 setInputData();
4037
4038 fitter.fit();
4039
4040 final var result = fitter.getA();
4041
4042 final var bx = result[0];
4043 final var by = result[1];
4044 final var bz = result[2];
4045
4046 final var m11 = result[3];
4047 final var m21 = result[4];
4048 final var m31 = result[5];
4049
4050 final var m12 = result[6];
4051 final var m22 = result[7];
4052 final var m32 = result[8];
4053
4054 final var m13 = result[9];
4055 final var m23 = result[10];
4056 final var m33 = result[11];
4057
4058 final var g11 = result[12];
4059 final var g21 = result[13];
4060 final var g31 = result[14];
4061
4062 final var g12 = result[15];
4063 final var g22 = result[16];
4064 final var g32 = result[17];
4065
4066 final var g13 = result[18];
4067 final var g23 = result[19];
4068 final var g33 = result[20];
4069
4070 final var b = new Matrix(BodyKinematics.COMPONENTS, 1);
4071 b.setElementAtIndex(0, bx);
4072 b.setElementAtIndex(1, by);
4073 b.setElementAtIndex(2, bz);
4074
4075 final var m = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
4076 m.setElementAtIndex(0, m11);
4077 m.setElementAtIndex(1, m21);
4078 m.setElementAtIndex(2, m31);
4079
4080 m.setElementAtIndex(3, m12);
4081 m.setElementAtIndex(4, m22);
4082 m.setElementAtIndex(5, m32);
4083
4084 m.setElementAtIndex(6, m13);
4085 m.setElementAtIndex(7, m23);
4086 m.setElementAtIndex(8, m33);
4087
4088 final var g = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
4089 g.setElementAtIndex(0, g11);
4090 g.setElementAtIndex(1, g21);
4091 g.setElementAtIndex(2, g31);
4092
4093 g.setElementAtIndex(3, g12);
4094 g.setElementAtIndex(4, g22);
4095 g.setElementAtIndex(5, g32);
4096
4097 g.setElementAtIndex(6, g13);
4098 g.setElementAtIndex(7, g23);
4099 g.setElementAtIndex(8, g33);
4100
4101 setResult(m, b, g);
4102
4103 // at this point covariance is expressed in terms of b, M and G, and must
4104 // be expressed in terms of bg, Mg and Gg.
4105 // We know that:
4106 // bg = M * b
4107 // Mg = M - I
4108 // Gg = M * G
4109
4110 // b = [bx]
4111 // [by]
4112 // [bz]
4113
4114 // M = [m11 m12 m13]
4115 // [m21 m22 m23]
4116 // [m31 m32 m33]
4117
4118 // G = [g11 g12 g13]
4119 // [g21 g22 g23]
4120 // [g31 g32 g33]
4121
4122 // bg = [m11 m12 m13][bx] = [m11 * bx + m12 * by + m13 * bz] = [bgx]
4123 // [m21 m22 m23][by] [m21 * bx + m22 * by + m23 * bz] [bgy]
4124 // [m31 m32 m33][bz] [m31 * bx + m32 * by + m33 * bz] [bgz]
4125
4126 // Mg = [sx mxy mxz] = [m11 - 1 m12 m13 ]
4127 // [myx sy myz] [m21 m22 - 1 m23 ]
4128 // [mzx mzy sz ] [m31 m32 m33 - 1 ]
4129
4130 // Gg = [gg11 gg12 gg13] = [m11 m12 m13][g11 g12 g13]
4131 // [gg21 gg22 gg23] [m21 m22 m23][g21 g22 g23]
4132 // [gg31 gg32 gg33] [m31 m32 m33][g31 g32 g33]
4133
4134 // Defining the linear application:
4135 // F(b, M, G) = F(bx, by, bz, m11, m21, m31, m12, m22, m32, m13, m23, m33, g11, g21, g31, g12, g22, g32, g13, g23, g33)
4136 // as:
4137 // [bgx] = [m11 * bx + m12 * by + m13 * bz]
4138 // [bgy] [m21 * bx + m22 * by + m23 * bz]
4139 // [bgz] [m31 * bx + m32 * by + m33 * bz]
4140 // [sx] [m11 - 1]
4141 // [sy] [m22 - 1]
4142 // [sz] [m33 - 1]
4143 // [mxy] [m12]
4144 // [mxz] [m13]
4145 // [myx] [m21]
4146 // [myz] [m23]
4147 // [mzx] [m31]
4148 // [mzy] [m32]
4149 // [gg11] [m11 * g11 + m12 * g21 + m13 * g31]
4150 // [gg21] [m21 * g11 + m22 * g21 + m23 * g31]
4151 // [gg31] [m31 * g11 + m32 * g21 + m33 * g31]
4152 // [gg12] [m11 * g12 + m12 * g22 + m13 * g32]
4153 // [gg22] [m21 * g12 + m22 * g22 + m23 * g32]
4154 // [gg32] [m31 * g12 + m32 * g22 + m33 * g32]
4155 // [gg13] [m11 * g13 + m12 * g23 + m13 * g33]
4156 // [gg23] [m21 * g13 + m22 * g23 + m23 * g33]
4157 // [gg33] [m31 * g13 + m32 * g23 + m33 * g33]
4158
4159 // Then the Jacobian of F(b, M, G) is:
4160 // J = [m11 m12 m13 bx 0 0 by 0 0 bz 0 0 0 0 0 0 0 0 0 0 0 ]
4161 // [m21 m22 m23 0 bx 0 0 by 0 0 bz 0 0 0 0 0 0 0 0 0 0 ]
4162 // [m31 m32 m33 0 0 bx 0 0 by 0 0 bz 0 0 0 0 0 0 0 0 0 ]
4163 // [0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ]
4164 // [0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 ]
4165 // [0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 ]
4166 // [0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ]
4167 // [0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 ]
4168 // [0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ]
4169 // [0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 ]
4170 // [0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ]
4171 // [0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 ]
4172 // [0 0 0 g11 0 0 g21 0 0 g31 0 0 m11 m12 m13 0 0 0 0 0 0 ]
4173 // [0 0 0 0 g11 0 0 g21 0 0 g31 0 m21 m22 m23 0 0 0 0 0 0 ]
4174 // [0 0 0 0 0 g11 0 0 g21 0 0 g31 m31 m32 m33 0 0 0 0 0 0 ]
4175 // [0 0 0 g12 0 0 g22 0 0 g32 0 0 0 0 0 m11 m12 m13 0 0 0 ]
4176 // [0 0 0 0 g12 0 0 g22 0 0 g32 0 0 0 0 m21 m22 m23 0 0 0 ]
4177 // [0 0 0 0 0 g12 0 0 g22 0 0 g32 0 0 0 m31 m32 m33 0 0 0 ]
4178 // [0 0 0 g13 0 0 g23 0 0 g33 0 0 0 0 0 0 0 0 m11 m12 m13]
4179 // [0 0 0 0 g13 0 0 g23 0 0 g33 0 0 0 0 0 0 0 m21 m22 m23]
4180 // [0 0 0 0 0 g13 0 0 g23 0 0 g33 0 0 0 0 0 0 m31 m32 m33]
4181
4182 // We know that the propagated covariance is J * Cov * J', hence:
4183 final var jacobian = new Matrix(GENERAL_UNKNOWNS_AND_CROSS_BIASES, GENERAL_UNKNOWNS_AND_CROSS_BIASES);
4184
4185 jacobian.setElementAt(0, 0, m11);
4186 jacobian.setElementAt(0, 1, m12);
4187 jacobian.setElementAt(0, 2, m13);
4188 jacobian.setElementAt(0, 3, bx);
4189 jacobian.setElementAt(0, 6, by);
4190 jacobian.setElementAt(0, 9, bz);
4191
4192 jacobian.setElementAt(1, 0, m21);
4193 jacobian.setElementAt(1, 1, m22);
4194 jacobian.setElementAt(1, 2, m23);
4195 jacobian.setElementAt(1, 4, bx);
4196 jacobian.setElementAt(1, 7, by);
4197 jacobian.setElementAt(1, 10, bz);
4198
4199 jacobian.setElementAt(2, 0, m31);
4200 jacobian.setElementAt(2, 1, m32);
4201 jacobian.setElementAt(2, 2, m33);
4202 jacobian.setElementAt(2, 5, bx);
4203 jacobian.setElementAt(2, 8, by);
4204 jacobian.setElementAt(2, 11, bz);
4205
4206 jacobian.setElementAt(3, 3, 1.0);
4207 jacobian.setElementAt(4, 7, 1.0);
4208 jacobian.setElementAt(5, 11, 1.0);
4209
4210 jacobian.setElementAt(6, 6, 1.0);
4211 jacobian.setElementAt(7, 9, 1.0);
4212 jacobian.setElementAt(8, 4, 1.0);
4213
4214 jacobian.setElementAt(9, 10, 1.0);
4215 jacobian.setElementAt(10, 5, 1.0);
4216 jacobian.setElementAt(11, 8, 1.9);
4217
4218 jacobian.setElementAt(12, 3, g11);
4219 jacobian.setElementAt(12, 6, g21);
4220 jacobian.setElementAt(12, 9, g31);
4221 jacobian.setElementAt(12, 12, m11);
4222 jacobian.setElementAt(12, 13, m12);
4223 jacobian.setElementAt(12, 14, m13);
4224
4225 jacobian.setElementAt(13, 4, g11);
4226 jacobian.setElementAt(13, 7, g21);
4227 jacobian.setElementAt(13, 10, g31);
4228 jacobian.setElementAt(13, 12, m21);
4229 jacobian.setElementAt(13, 13, m22);
4230 jacobian.setElementAt(13, 14, m23);
4231
4232 jacobian.setElementAt(14, 5, g11);
4233 jacobian.setElementAt(14, 8, g21);
4234 jacobian.setElementAt(14, 11, g31);
4235 jacobian.setElementAt(14, 12, m31);
4236 jacobian.setElementAt(14, 13, m32);
4237 jacobian.setElementAt(14, 14, m33);
4238
4239 jacobian.setElementAt(15, 3, g12);
4240 jacobian.setElementAt(15, 6, g22);
4241 jacobian.setElementAt(15, 9, g32);
4242 jacobian.setElementAt(15, 15, m11);
4243 jacobian.setElementAt(15, 16, m12);
4244 jacobian.setElementAt(15, 17, m13);
4245
4246 jacobian.setElementAt(16, 4, g12);
4247 jacobian.setElementAt(16, 7, g22);
4248 jacobian.setElementAt(16, 10, g32);
4249 jacobian.setElementAt(16, 15, m21);
4250 jacobian.setElementAt(16, 16, m22);
4251 jacobian.setElementAt(16, 17, m23);
4252
4253 jacobian.setElementAt(17, 5, g12);
4254 jacobian.setElementAt(17, 8, g22);
4255 jacobian.setElementAt(17, 11, g32);
4256 jacobian.setElementAt(17, 15, m31);
4257 jacobian.setElementAt(17, 16, m32);
4258 jacobian.setElementAt(17, 17, m33);
4259
4260 jacobian.setElementAt(18, 3, g13);
4261 jacobian.setElementAt(18, 6, g23);
4262 jacobian.setElementAt(18, 9, g33);
4263 jacobian.setElementAt(18, 18, m11);
4264 jacobian.setElementAt(18, 19, m12);
4265 jacobian.setElementAt(18, 20, m13);
4266
4267 jacobian.setElementAt(19, 4, g13);
4268 jacobian.setElementAt(19, 7, g23);
4269 jacobian.setElementAt(19, 10, g33);
4270 jacobian.setElementAt(19, 18, m21);
4271 jacobian.setElementAt(19, 19, m22);
4272 jacobian.setElementAt(19, 20, m23);
4273
4274 jacobian.setElementAt(20, 5, g13);
4275 jacobian.setElementAt(20, 8, g23);
4276 jacobian.setElementAt(20, 11, g33);
4277 jacobian.setElementAt(20, 18, m31);
4278 jacobian.setElementAt(20, 19, m32);
4279 jacobian.setElementAt(20, 20, m33);
4280
4281 final var jacobianTrans = jacobian.transposeAndReturnNew();
4282 jacobian.multiply(estimatedCovariance);
4283 jacobian.multiply(jacobianTrans);
4284 estimatedCovariance = jacobian;
4285 }
4286
4287 /**
4288 * Internal method to perform calibration when common z-axis is assumed
4289 * for both the accelerometer and gyroscope and G-dependent cross biases
4290 * are ignored.
4291 *
4292 * @throws AlgebraException if accelerometer parameters prevent fixing
4293 * measured accelerometer values due to numerical instabilities.
4294 * @throws FittingException if no convergence to solution is found.
4295 * @throws com.irurueta.numerical.NotReadyException if fitter is not ready.
4296 */
4297 private void calibrateCommonAxis() throws AlgebraException, FittingException,
4298 com.irurueta.numerical.NotReadyException {
4299 // The gyroscope model is
4300 // Ωmeas = bg + (I + Mg) * Ωtrue + Gg * ftrue + w
4301
4302 // Ideally a least squares solution tries to minimize noise component, so:
4303 // Ωmeas = bg + (I + Mg) * Ωtrue + Gg * ftrue
4304
4305 // Where Mg is upper triangular
4306
4307 // Since G-dependent cross biases are ignored, we can assume that Gg = 0
4308
4309 // Hence:
4310 // Ωmeas = bg + (I + Mg) * Ωtrue
4311
4312 // For convergence purposes of the Levenberg-Marquardt algorithm,
4313 // the gyroscope model can be better expressed as:
4314 // Ωmeas = T*K*(Ωtrue + b)
4315 // Ωmeas = M*(Ωtrue + b)
4316 // Ωmeas = M*Ωtrue + M*b
4317
4318 // Hence
4319 // Ωmeas - M*b = M*Ωtrue
4320
4321 // M^-1 * (Ωmeas - M*b) = Ωtrue
4322
4323 // where:
4324 // M = I + Mg
4325 // bg = M*b = (I + Mg)*b --> b = M^-1*bg
4326
4327 // Notice that M is upper diagonal because Mg is upper diagonal
4328 // when common axis is assumed
4329
4330 final var gradientEstimator = new GradientEstimator(params -> evaluateCommonAxis(i, params));
4331
4332 final var initialM = Matrix.identity(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
4333 initialM.add(getInitialMg());
4334
4335 // Force initial M to be upper diagonal
4336 initialM.setElementAt(1, 0, 0.0);
4337 initialM.setElementAt(2, 0, 0.0);
4338 initialM.setElementAt(2, 1, 0.0);
4339
4340 final var invInitialM = Utils.inverse(initialM);
4341 final var initialBg = getInitialBiasAsMatrix();
4342 final var initialB = invInitialM.multiplyAndReturnNew(initialBg);
4343
4344 fitter.setFunctionEvaluator(new LevenbergMarquardtMultiDimensionFunctionEvaluator() {
4345 @Override
4346 public int getNumberOfDimensions() {
4347 // Before and after normalized gravity versors
4348 return 2 * BodyKinematics.COMPONENTS;
4349 }
4350
4351 @Override
4352 public double[] createInitialParametersArray() {
4353 final var initial = new double[COMMON_Z_AXIS_UNKNOWNS];
4354
4355 // biases b
4356 for (var i = 0; i < BodyKinematics.COMPONENTS; i++) {
4357 initial[i] = initialB.getElementAtIndex(i);
4358 }
4359
4360 // upper diagonal cross coupling errors M
4361 var k = BodyKinematics.COMPONENTS;
4362 for (var j = 0; j < BodyKinematics.COMPONENTS; j++) {
4363 for (var i = 0; i < BodyKinematics.COMPONENTS; i++) {
4364 if (i <= j) {
4365 initial[k] = initialM.getElementAt(i, j);
4366 k++;
4367 }
4368 }
4369 }
4370
4371 return initial;
4372 }
4373
4374 @Override
4375 public double evaluate(
4376 final int i, final double[] point, final double[] params, final double[] derivatives)
4377 throws EvaluationException {
4378 EasyGyroscopeCalibrator.this.i = i;
4379
4380 // point contains fixed gravity versor values for current
4381 // sequence
4382 EasyGyroscopeCalibrator.this.point = point;
4383
4384 gradientEstimator.gradient(params, derivatives);
4385
4386 return evaluateCommonAxis(i, params);
4387 }
4388 });
4389
4390 setInputData();
4391
4392 fitter.fit();
4393
4394 final var result = fitter.getA();
4395
4396 final var bx = result[0];
4397 final var by = result[1];
4398 final var bz = result[2];
4399
4400 final var m11 = result[3];
4401
4402 final var m12 = result[4];
4403 final var m22 = result[5];
4404
4405 final var m13 = result[6];
4406 final var m23 = result[7];
4407 final var m33 = result[8];
4408
4409 final var b = new Matrix(BodyKinematics.COMPONENTS, 1);
4410 b.setElementAtIndex(0, bx);
4411 b.setElementAtIndex(1, by);
4412 b.setElementAtIndex(2, bz);
4413
4414 final var m = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
4415 m.setElementAtIndex(0, m11);
4416 m.setElementAtIndex(1, 0.0);
4417 m.setElementAtIndex(2, 0.0);
4418
4419 m.setElementAtIndex(3, m12);
4420 m.setElementAtIndex(4, m22);
4421 m.setElementAtIndex(5, 0.0);
4422
4423 m.setElementAtIndex(6, m13);
4424 m.setElementAtIndex(7, m23);
4425 m.setElementAtIndex(8, m33);
4426
4427 setResult(m, b);
4428
4429 // at this point covariance is expressed in terms of b, M and G, and must
4430 // be expressed in terms of bg, Mg and Gg.
4431 // We know that:
4432 // bg = M * b
4433 // Mg = M - I
4434 // Gg = M * G = 0
4435
4436 // b = [bx]
4437 // [by]
4438 // [bz]
4439
4440 // M = [m11 m12 m13]
4441 // [0 m22 m23]
4442 // [0 0 m33]
4443
4444 // G = [g11 g12 g13] = 0
4445 // [g21 g22 g23]
4446 // [g31 g32 g33]
4447
4448 // bg = [m11 m12 m13][bx] = [m11 * bx + m12 * by + m13 * bz] = [bgx]
4449 // [0 m22 m23][by] [ m22 * by + m23 * bz] [bgy]
4450 // [0 0 m33][bz] [ m33 * bz] [bgz]
4451
4452 // Mg = [sx mxy mxz] = [m11 - 1 m12 m13 ]
4453 // [myx sy myz] [0 m22 - 1 m23 ]
4454 // [mzx mzy sz ] [0 0 m33 - 1 ]
4455
4456 // Gg = [gg11 gg12 gg13] = 0
4457 // [gg21 gg22 gg23]
4458 // [gg31 gg32 gg33]
4459
4460 // Defining the linear application:
4461 // F(b, M) = F(bx, by, bz, m11, m12, m22, m13, m23, m33)
4462 // as:
4463 // [bgx] = [m11 * bx + m12 * by + m13 * bz]
4464 // [bgy] [ m22 * by + m23 * bz]
4465 // [bgz] [ m33 * bz]
4466 // [sx] [m11 - 1]
4467 // [sy] [m22 - 1]
4468 // [sz] [m33 - 1]
4469 // [mxy] [m12]
4470 // [mxz] [m13]
4471 // [myx] [0]
4472 // [myz] [m23]
4473 // [mzx] [0]
4474 // [mzy] [0]
4475 // [gg11] [0]
4476 // [gg21] [0]
4477 // [gg31] [0]
4478 // [gg12] [0]
4479 // [gg22] [0]
4480 // [gg32] [0]
4481 // [gg13] [0]
4482 // [gg23] [0]
4483 // [gg33] [0]
4484
4485 // Then the Jacobian of F(b, M) is:
4486 // J = [m11 m12 m13 bx by 0 bz 0 0 ]
4487 // [0 m22 m23 0 0 by 0 bz 0 ]
4488 // [0 0 m33 0 0 0 0 0 bz ]
4489 // [0 0 0 1 0 0 0 0 0 ]
4490 // [0 0 0 0 0 1 0 0 0 ]
4491 // [0 0 0 0 0 0 0 0 1 ]
4492 // [0 0 0 0 1 0 0 0 0 ]
4493 // [0 0 0 0 0 0 1 0 0 ]
4494 // [0 0 0 0 0 0 0 0 0 ]
4495 // [0 0 0 0 0 0 0 1 0 ]
4496 // [0 0 0 0 0 0 0 0 0 ]
4497 // [0 0 0 0 0 0 0 0 0 ]
4498 // [0 0 0 0 0 0 0 0 0 ]
4499 // [0 0 0 0 0 0 0 0 0 ]
4500 // [0 0 0 0 0 0 0 0 0 ]
4501 // [0 0 0 0 0 0 0 0 0 ]
4502 // [0 0 0 0 0 0 0 0 0 ]
4503 // [0 0 0 0 0 0 0 0 0 ]
4504 // [0 0 0 0 0 0 0 0 0 ]
4505 // [0 0 0 0 0 0 0 0 0 ]
4506 // [0 0 0 0 0 0 0 0 0 ]
4507
4508 // We know that the propagated covariance is J * Cov * J', hence:
4509 final var jacobian = new Matrix(GENERAL_UNKNOWNS_AND_CROSS_BIASES, COMMON_Z_AXIS_UNKNOWNS);
4510
4511 jacobian.setElementAt(0, 0, m11);
4512 jacobian.setElementAt(0, 1, m12);
4513 jacobian.setElementAt(0, 2, m13);
4514 jacobian.setElementAt(0, 3, bx);
4515 jacobian.setElementAt(0, 4, by);
4516 jacobian.setElementAt(0, 6, bz);
4517
4518 jacobian.setElementAt(1, 1, m22);
4519 jacobian.setElementAt(1, 2, m23);
4520 jacobian.setElementAt(1, 5, by);
4521 jacobian.setElementAt(1, 7, bz);
4522
4523 jacobian.setElementAt(2, 2, m33);
4524 jacobian.setElementAt(2, 8, bz);
4525
4526 jacobian.setElementAt(3, 3, 1.0);
4527 jacobian.setElementAt(4, 5, 1.0);
4528 jacobian.setElementAt(5, 8, 1.0);
4529
4530 jacobian.setElementAt(6, 4, 1.0);
4531 jacobian.setElementAt(7, 6, 1.0);
4532
4533 jacobian.setElementAt(9, 7, 1.0);
4534
4535 final var jacobianTrans = jacobian.transposeAndReturnNew();
4536 jacobian.multiply(estimatedCovariance);
4537 jacobian.multiply(jacobianTrans);
4538 estimatedCovariance = jacobian;
4539 }
4540
4541 /**
4542 * Internal method to perform general calibration when G-dependant cross
4543 * biases are ignored.
4544 *
4545 * @throws AlgebraException if accelerometer parameters prevent fixing
4546 * measured accelerometer values due to numerical instabilities.
4547 * @throws FittingException if no convergence to solution is found.
4548 * @throws com.irurueta.numerical.NotReadyException if fitter is not ready.
4549 */
4550 private void calibrateGeneral() throws AlgebraException, FittingException,
4551 com.irurueta.numerical.NotReadyException {
4552 // The gyroscope model is
4553 // Ωmeas = bg + (I + Mg) * Ωtrue + Gg * ftrue + w
4554
4555 // Ideally a least squares solution tries to minimize noise component, so:
4556 // Ωmeas = bg + (I + Mg) * Ωtrue + Gg * ftrue
4557
4558 // Since G-dependent cross biases are ignored, we can assume that Gg = 0
4559
4560 // Hence:
4561 // Ωmeas = bg + (I + Mg) * Ωtrue
4562
4563 // For convergence purposes of the Levenberg-Marquardt algorithm, the
4564 // gyroscope model can be better expressed as:
4565 // Ωmeas = T*K*(Ωtrue + b)
4566 // Ωmeas = M*(Ωtrue + b)
4567 // Ωmeas = M*Ωtrue + M*b
4568
4569 // where:
4570 // M = I + Mg
4571 // bg = M*b = (I + Mg)*b --> b = M^-1*bg
4572
4573 // Hence
4574 // Ωmeas - M*b = M*Ωtrue
4575
4576 // M^-1 * (Ωmeas - M*b) = Ωtrue
4577
4578 final var gradientEstimator = new GradientEstimator(params -> evaluateGeneral(i, params));
4579
4580 final var initialM = Matrix.identity(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
4581 initialM.add(getInitialMg());
4582
4583 final var invInitialM = Utils.inverse(initialM);
4584 final var initialBg = getInitialBiasAsMatrix();
4585 final var initialB = invInitialM.multiplyAndReturnNew(initialBg);
4586
4587 fitter.setFunctionEvaluator(new LevenbergMarquardtMultiDimensionFunctionEvaluator() {
4588 @Override
4589 public int getNumberOfDimensions() {
4590 // Before and after normalized gravity versors
4591 return 2 * BodyKinematics.COMPONENTS;
4592 }
4593
4594 @Override
4595 public double[] createInitialParametersArray() {
4596 final var initial = new double[GENERAL_UNKNOWNS];
4597
4598 // biases b
4599 for (var i = 0; i < BodyKinematics.COMPONENTS; i++) {
4600 initial[i] = initialB.getElementAtIndex(i);
4601 }
4602
4603 // cross coupling errors M
4604 final var num = BodyKinematics.COMPONENTS * BodyKinematics.COMPONENTS;
4605 for (int i = 0, j = BodyKinematics.COMPONENTS; i < num; i++, j++) {
4606 initial[j] = initialM.getElementAtIndex(i);
4607 }
4608
4609 return initial;
4610 }
4611
4612 @Override
4613 public double evaluate(
4614 final int i, final double[] point, final double[] params, final double[] derivatives)
4615 throws EvaluationException {
4616 EasyGyroscopeCalibrator.this.i = i;
4617
4618 // point contains fixed gravity versor values for current
4619 // sequence
4620 EasyGyroscopeCalibrator.this.point = point;
4621
4622 gradientEstimator.gradient(params, derivatives);
4623
4624 return evaluateGeneral(i, params);
4625 }
4626 });
4627
4628 setInputData();
4629
4630 fitter.fit();
4631
4632 final var result = fitter.getA();
4633
4634 final var bx = result[0];
4635 final var by = result[1];
4636 final var bz = result[2];
4637
4638 final var m11 = result[3];
4639 final var m21 = result[4];
4640 final var m31 = result[5];
4641
4642 final var m12 = result[6];
4643 final var m22 = result[7];
4644 final var m32 = result[8];
4645
4646 final var m13 = result[9];
4647 final var m23 = result[10];
4648 final var m33 = result[11];
4649
4650 final var b = new Matrix(BodyKinematics.COMPONENTS, 1);
4651 b.setElementAtIndex(0, bx);
4652 b.setElementAtIndex(1, by);
4653 b.setElementAtIndex(2, bz);
4654
4655 final var m = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
4656 m.setElementAtIndex(0, m11);
4657 m.setElementAtIndex(1, m21);
4658 m.setElementAtIndex(2, m31);
4659
4660 m.setElementAtIndex(3, m12);
4661 m.setElementAtIndex(4, m22);
4662 m.setElementAtIndex(5, m32);
4663
4664 m.setElementAtIndex(6, m13);
4665 m.setElementAtIndex(7, m23);
4666 m.setElementAtIndex(8, m33);
4667
4668 setResult(m, b);
4669
4670 // at this point covariance is expressed in terms of b, M and G, and must
4671 // be expressed in terms of bg, Mg and Gg.
4672 // We know that:
4673 // bg = M * b
4674 // Mg = M - I
4675 // Gg = M * G = 0
4676
4677 // b = [bx]
4678 // [by]
4679 // [bz]
4680
4681 // M = [m11 m12 m13]
4682 // [m21 m22 m23]
4683 // [m31 m32 m33]
4684
4685 // G = [g11 g12 g13] = 0
4686 // [g21 g22 g23]
4687 // [g31 g32 g33]
4688
4689 // bg = [m11 m12 m13][bx] = [m11 * bx + m12 * by + m13 * bz] = [bgx]
4690 // [m21 m22 m23][by] [m21 * bx + m22 * by + m23 * bz] [bgy]
4691 // [m31 m32 m33][bz] [m31 * bx + m32 * by + m33 * bz] [bgz]
4692
4693 // Mg = [sx mxy mxz] = [m11 - 1 m12 m13 ]
4694 // [myx sy myz] [m21 m22 - 1 m23 ]
4695 // [mzx mzy sz ] [m31 m32 m33 - 1 ]
4696
4697 // Gg = [gg11 gg12 gg13] = 0
4698 // [gg21 gg22 gg23]
4699 // [gg31 gg32 gg33]
4700
4701 // Defining the linear application:
4702 // F(b, M) = F(bx, by, bz, m11, m21, m31, m12, m22, m32, m13, m23, m33)
4703 // as:
4704 // [bgx] = [m11 * bx + m12 * by + m13 * bz]
4705 // [bgy] [m21 * bx + m22 * by + m23 * bz]
4706 // [bgz] [m31 * bx + m32 * by + m33 * bz]
4707 // [sx] [m11 - 1]
4708 // [sy] [m22 - 1]
4709 // [sz] [m33 - 1]
4710 // [mxy] [m12]
4711 // [mxz] [m13]
4712 // [myx] [m21]
4713 // [myz] [m23]
4714 // [mzx] [m31]
4715 // [mzy] [m32]
4716 // [gg11] [0]
4717 // [gg21] [0]
4718 // [gg31] [0]
4719 // [gg12] [0]
4720 // [gg22] [0]
4721 // [gg32] [0]
4722 // [gg13] [0]
4723 // [gg23] [0]
4724 // [gg33] [0]
4725
4726 // Then the Jacobian of F(b, M) is:
4727 // J = [m11 m12 m13 bx 0 0 by 0 0 bz 0 0 ]
4728 // [m21 m22 m23 0 bx 0 0 by 0 0 bz 0 ]
4729 // [m31 m32 m33 0 0 bx 0 0 by 0 0 bz]
4730 // [0 0 0 1 0 0 0 0 0 0 0 0 ]
4731 // [0 0 0 0 0 0 0 1 0 0 0 0 ]
4732 // [0 0 0 0 0 0 0 0 0 0 0 1 ]
4733 // [0 0 0 0 0 0 1 0 0 0 0 0 ]
4734 // [0 0 0 0 0 0 0 0 0 1 0 0 ]
4735 // [0 0 0 0 1 0 0 0 0 0 0 0 ]
4736 // [0 0 0 0 0 0 0 0 0 0 1 0 ]
4737 // [0 0 0 0 0 1 0 0 0 0 0 0 ]
4738 // [0 0 0 0 0 0 0 0 1 0 0 0 ]
4739 // [0 0 0 0 0 0 0 0 0 0 0 0 ]
4740 // [0 0 0 0 0 0 0 0 0 0 0 0 ]
4741 // [0 0 0 0 0 0 0 0 0 0 0 0 ]
4742 // [0 0 0 0 0 0 0 0 0 0 0 0 ]
4743 // [0 0 0 0 0 0 0 0 0 0 0 0 ]
4744 // [0 0 0 0 0 0 0 0 0 0 0 0 ]
4745 // [0 0 0 0 0 0 0 0 0 0 0 0 ]
4746 // [0 0 0 0 0 0 0 0 0 0 0 0 ]
4747 // [0 0 0 0 0 0 0 0 0 0 0 0 ]
4748
4749 // We know that the propagated covariance is J * Cov * J', hence:
4750 final var jacobian = new Matrix(GENERAL_UNKNOWNS_AND_CROSS_BIASES, GENERAL_UNKNOWNS);
4751
4752 jacobian.setElementAt(0, 0, m11);
4753 jacobian.setElementAt(0, 1, m12);
4754 jacobian.setElementAt(0, 2, m13);
4755 jacobian.setElementAt(0, 3, bx);
4756 jacobian.setElementAt(0, 6, by);
4757 jacobian.setElementAt(0, 9, bz);
4758
4759 jacobian.setElementAt(1, 0, m21);
4760 jacobian.setElementAt(1, 1, m22);
4761 jacobian.setElementAt(1, 2, m23);
4762 jacobian.setElementAt(1, 4, bx);
4763 jacobian.setElementAt(1, 7, by);
4764 jacobian.setElementAt(1, 10, bz);
4765
4766 jacobian.setElementAt(2, 0, m31);
4767 jacobian.setElementAt(2, 1, m32);
4768 jacobian.setElementAt(2, 2, m33);
4769 jacobian.setElementAt(2, 5, bx);
4770 jacobian.setElementAt(2, 8, by);
4771 jacobian.setElementAt(2, 11, bz);
4772
4773 jacobian.setElementAt(3, 3, 1.0);
4774 jacobian.setElementAt(4, 7, 1.0);
4775 jacobian.setElementAt(5, 11, 1.0);
4776
4777 jacobian.setElementAt(6, 6, 1.0);
4778 jacobian.setElementAt(7, 9, 1.0);
4779 jacobian.setElementAt(8, 4, 1.0);
4780
4781 jacobian.setElementAt(9, 10, 1.0);
4782 jacobian.setElementAt(10, 5, 1.0);
4783 jacobian.setElementAt(11, 8, 1.9);
4784
4785 final var jacobianTrans = jacobian.transposeAndReturnNew();
4786 jacobian.multiply(estimatedCovariance);
4787 jacobian.multiply(jacobianTrans);
4788 estimatedCovariance = jacobian;
4789 }
4790
4791 /**
4792 * Sets input data into Levenberg-Marquardt fitter.
4793 *
4794 * @throws AlgebraException if there are numerical instabilities.
4795 */
4796 private void setInputData() throws AlgebraException {
4797
4798 final var ba = getAccelerometerBiasAsMatrix();
4799 final var ma = getAccelerometerMa();
4800
4801 final var measuredBeforeF = new double[BodyKinematics.COMPONENTS];
4802 final var fixedBeforeF = new double[BodyKinematics.COMPONENTS];
4803
4804 final var measuredAfterF = new double[BodyKinematics.COMPONENTS];
4805 final var fixedAfterF = new double[BodyKinematics.COMPONENTS];
4806
4807 final var numSequences = sequences.size();
4808 final var x = new Matrix(numSequences, 2 * BodyKinematics.COMPONENTS);
4809 final var y = new double[numSequences];
4810 final var standardDeviations = new double[numSequences];
4811
4812 // make a copy of input sequences that will be used to update
4813 // kinematics measurements with fixed values for memory efficiency
4814
4815 fixedSequences = new ArrayList<>();
4816 for (final var sequence : sequences) {
4817 fixedSequences.add(new BodyKinematicsSequence<>(sequence));
4818 }
4819
4820 accelerationFixer.setBias(ba);
4821 accelerationFixer.setCrossCouplingErrors(ma);
4822
4823 var i = 0;
4824 for (final var sequence : sequences) {
4825 // sequence mean accelerometer samples of previous static
4826 // period will need to be fixed using accelerometer calibration
4827 // parameters
4828 measuredBeforeF[0] = sequence.getBeforeMeanFx();
4829 measuredBeforeF[1] = sequence.getBeforeMeanFy();
4830 measuredBeforeF[2] = sequence.getBeforeMeanFz();
4831 accelerationFixer.fix(measuredBeforeF, fixedBeforeF);
4832
4833 measuredAfterF[0] = sequence.getAfterMeanFx();
4834 measuredAfterF[1] = sequence.getAfterMeanFy();
4835 measuredAfterF[2] = sequence.getAfterMeanFz();
4836 accelerationFixer.fix(measuredAfterF, fixedAfterF);
4837
4838 // because we are only interested in gravity direction, we
4839 // normalize these vectors, so that gravity becomes independent
4840 // of current Earth position.
4841 ArrayUtils.normalize(fixedBeforeF);
4842 ArrayUtils.normalize(fixedAfterF);
4843
4844 x.setSubmatrix(i, 0, i, 2, fixedBeforeF);
4845 x.setSubmatrix(i, 3, i, 5, fixedAfterF);
4846
4847 y[i] = 0.0;
4848
4849 standardDeviations[i] = computeAverageAngularRateStandardDeviation(sequence);
4850 i++;
4851 }
4852
4853 fitter.setInputData(x, y, standardDeviations);
4854 }
4855
4856 /**
4857 * Computes average angular rate standard deviation for measurements
4858 * in provided sequence.
4859 *
4860 * @param sequence a sequence.
4861 * @return average angular rate standard deviation expressed in radians
4862 * per second (rad/s).
4863 */
4864 private static double computeAverageAngularRateStandardDeviation(
4865 final BodyKinematicsSequence<StandardDeviationTimedBodyKinematics> sequence) {
4866 final var items = sequence.getSortedItems();
4867 final var size = items.size();
4868
4869 var result = 0.0;
4870 for (final var item : items) {
4871 result += item.getAngularRateStandardDeviation() / size;
4872 }
4873
4874 return result;
4875 }
4876
4877 /**
4878 * Converts acceleration instance to meters per squared second.
4879 *
4880 * @param acceleration acceleration instance to be converted.
4881 * @return converted value.
4882 */
4883 private static double convertAcceleration(final Acceleration acceleration) {
4884 return AccelerationConverter.convert(acceleration.getValue().doubleValue(), acceleration.getUnit(),
4885 AccelerationUnit.METERS_PER_SQUARED_SECOND);
4886 }
4887
4888 /**
4889 * Converts angular speed value and unit to radians per second.
4890 *
4891 * @param value angular speed value.
4892 * @param unit unit of angular speed value.
4893 * @return converted value.
4894 */
4895 private static double convertAngularSpeed(final double value, final AngularSpeedUnit unit) {
4896 return AngularSpeedConverter.convert(value, unit, AngularSpeedUnit.RADIANS_PER_SECOND);
4897 }
4898
4899 /**
4900 * Converts angular speed instance to radians per second.
4901 *
4902 * @param angularSpeed angular speed instance to be converted.
4903 * @return converted value.
4904 */
4905 private static double convertAngularSpeed(final AngularSpeed angularSpeed) {
4906 return convertAngularSpeed(angularSpeed.getValue().doubleValue(), angularSpeed.getUnit());
4907 }
4908
4909 /**
4910 * Makes proper conversion of internal cross-coupling, bias and g-dependent
4911 * cross bias matrices.
4912 *
4913 * @param m internal scaling and cross-coupling matrix.
4914 * @param b internal bias matrix.
4915 * @param g internal g-dependent cross bias matrix.
4916 * @throws AlgebraException if a numerical instability occurs.
4917 */
4918 private void setResult(final Matrix m, final Matrix b, final Matrix g) throws AlgebraException {
4919 setResult(m, b);
4920
4921 // Gg = M*G
4922 m.multiply(g, estimatedGg);
4923 }
4924
4925 /**
4926 * Makes proper conversion of internal cross-coupling and bias matrices.
4927 *
4928 * @param m internal scaling and cross-coupling matrix.
4929 * @param b internal bias matrix.
4930 * @throws AlgebraException if a numerical instability occurs.
4931 */
4932 private void setResult(final Matrix m, final Matrix b) throws AlgebraException {
4933 // Because:
4934 // M = I + Mg
4935 // b = M^-1*bg
4936
4937 // Then:
4938 // Mg = M - I
4939 // bg = M*b
4940
4941 if (estimatedBiases == null) {
4942 estimatedBiases = new double[BodyKinematics.COMPONENTS];
4943 }
4944
4945 final var bg = m.multiplyAndReturnNew(b);
4946 bg.toArray(estimatedBiases);
4947
4948 if (estimatedMg == null) {
4949 estimatedMg = m;
4950 } else {
4951 estimatedMg.copyFrom(m);
4952 }
4953
4954 for (var i = 0; i < BodyKinematics.COMPONENTS; i++) {
4955 estimatedMg.setElementAt(i, i, estimatedMg.getElementAt(i, i) - 1.0);
4956 }
4957
4958 if (estimatedGg == null) {
4959 estimatedGg = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
4960 } else {
4961 estimatedGg.initialize(0.0);
4962 }
4963
4964 estimatedCovariance = fitter.getCovar();
4965 estimatedChiSq = fitter.getChisq();
4966 estimatedChiSqDegreesOfFreedom = fitter.getChisqDegreesOfFreedom();
4967 estimatedReducedChiSq = fitter.getReducedChisq();
4968 estimatedMse = fitter.getMse();
4969 try {
4970 estimatedP = fitter.getP();
4971 estimatedQ = fitter.getQ();
4972 } catch (final MaxIterationsExceededException ignore) {
4973 // if numerical instabilities arise, we assume worst case (no fit at all)
4974 // probability of finding a smaller chi square value is 1.0
4975 // quality of fit is 0.0
4976 estimatedP = 1.0;
4977 estimatedQ = 0.0;
4978 }
4979 }
4980
4981 /**
4982 * Computes gravity versor error at the end of a sequence using provided
4983 * parameters.
4984 * This method is internally executed during gradient estimation and
4985 * Levenberg-Marquardt fitting needed for calibration computation.
4986 *
4987 * @param i row position.
4988 * @param params array containing parameters for the general purpose case
4989 * when G-dependent cross biases are taken into account. Must
4990 * have length 18.
4991 * @return error between estimated and measured gravity versor.
4992 * @throws EvaluationException if there are numerical instabilities.
4993 */
4994 private double evaluateGeneralWithGDependentCrossBiases(final int i, final double[] params)
4995 throws EvaluationException {
4996
4997 final var bx = params[0];
4998 final var by = params[1];
4999 final var bz = params[2];
5000
5001 final var m11 = params[3];
5002 final var m21 = params[4];
5003 final var m31 = params[5];
5004
5005 final var m12 = params[6];
5006 final var m22 = params[7];
5007 final var m32 = params[8];
5008
5009 final var m13 = params[9];
5010 final var m23 = params[10];
5011 final var m33 = params[11];
5012
5013 final var g11 = params[12];
5014 final var g21 = params[13];
5015 final var g31 = params[14];
5016
5017 final var g12 = params[15];
5018 final var g22 = params[16];
5019 final var g32 = params[17];
5020
5021 final var g13 = params[18];
5022 final var g23 = params[19];
5023 final var g33 = params[20];
5024
5025 return evaluate(i, bx, by, bz, m11, m21, m31, m12, m22, m32, m13, m23, m33,
5026 g11, g21, g31, g12, g22, g32, g13, g23, g33);
5027 }
5028
5029 /**
5030 * Computes gravity versor error at the end of a sequence using provided
5031 * parameters.
5032 * This method is internally executed during gradient estimation and
5033 * Levenberg-Marquardt fitting needed for calibration computation.
5034 *
5035 * @param i row position.
5036 * @param params array containing parameters for the general purpose case
5037 * when G-dependent cross biases are taken into account. Must
5038 * have length 15.
5039 * @return error between estimated and measured gravity versor.
5040 * @throws EvaluationException if there are numerical instabilities.
5041 */
5042 private double evaluateCommonAxisWithGDependentCrossBiases(final int i, final double[] params)
5043 throws EvaluationException {
5044
5045 final var bx = params[0];
5046 final var by = params[1];
5047 final var bz = params[2];
5048
5049 final var m11 = params[3];
5050
5051 final var m12 = params[4];
5052 final var m22 = params[5];
5053
5054 final var m13 = params[6];
5055 final var m23 = params[7];
5056 final var m33 = params[8];
5057
5058 final var g11 = params[9];
5059 final var g21 = params[10];
5060 final var g31 = params[11];
5061
5062 final var g12 = params[12];
5063 final var g22 = params[13];
5064 final var g32 = params[14];
5065
5066 final var g13 = params[15];
5067 final var g23 = params[16];
5068 final var g33 = params[17];
5069
5070 return evaluate(i, bx, by, bz, m11, 0.0, 0.0, m12, m22, 0.0, m13, m23, m33,
5071 g11, g21, g31, g12, g22, g32, g13, g23, g33);
5072 }
5073
5074 /**
5075 * Computes gravity versor error at the end of a sequence using provided
5076 * parameters.
5077 * This method is internally executed during gradient estimation and
5078 * Levenberg-Marquardt fitting needed for calibration computation.
5079 *
5080 * @param i row position.
5081 * @param params array containing current parameters for the general purpose case
5082 * when G-dependent cross biases are ignored. Must have length 9.
5083 * @return error between estimated and measured gravity versor.
5084 * @throws EvaluationException if there are numerical instabilities.
5085 */
5086 private double evaluateGeneral(final int i, final double[] params) throws EvaluationException {
5087
5088 final var bx = params[0];
5089 final var by = params[1];
5090 final var bz = params[2];
5091
5092 final var m11 = params[3];
5093 final var m21 = params[4];
5094 final var m31 = params[5];
5095
5096 final var m12 = params[6];
5097 final var m22 = params[7];
5098 final var m32 = params[8];
5099
5100 final var m13 = params[9];
5101 final var m23 = params[10];
5102 final var m33 = params[11];
5103
5104 return evaluate(i, bx, by, bz, m11, m21, m31, m12, m22, m32, m13, m23, m33,
5105 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0);
5106 }
5107
5108 /**
5109 * Computes gravity versor error at the end of a sequence using provided
5110 * parameters.
5111 * This method is internally executed during gradient estimation and
5112 * Levenberg-Marquardt fitting needed for calibration computation.
5113 *
5114 * @param i row position.
5115 * @param params array containing current parameters for the common z-axis case
5116 * when G-dependent cross biases are ignored. Must have length 9.
5117 * @return error between estimated and measured gravity versor.
5118 * @throws EvaluationException if there are numerical instabilities.
5119 */
5120 private double evaluateCommonAxis(final int i, final double[] params) throws EvaluationException {
5121 final var bx = params[0];
5122 final var by = params[1];
5123 final var bz = params[2];
5124
5125 final var m11 = params[3];
5126
5127 final var m12 = params[4];
5128 final var m22 = params[5];
5129
5130 final var m13 = params[6];
5131 final var m23 = params[7];
5132 final var m33 = params[8];
5133
5134 return evaluate(i, bx, by, bz, m11, 0.0, 0.0, m12, m22, 0.0, m13, m23, m33,
5135 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0);
5136 }
5137
5138 /**
5139 * Computes gravity versor error at the end of a sequence using provided
5140 * parameters.
5141 * This method is internally executed during gradient estimation and
5142 * Levenberg-Marquardt fitting needed for calibration computation.
5143 *
5144 * @param i row position.
5145 * @param bx x coordinate of bias.
5146 * @param by y coordinate of bias.
5147 * @param bz z coordinate of bias.
5148 * @param m11 element 1,1 of cross-coupling error matrix.
5149 * @param m21 element 2,1 of cross-coupling error matrix.
5150 * @param m31 element 3,1 of cross-coupling error matrix.
5151 * @param m12 element 1,2 of cross-coupling error matrix.
5152 * @param m22 element 2,2 of cross-coupling error matrix.
5153 * @param m32 element 3,2 of cross-coupling error matrix.
5154 * @param m13 element 1,3 of cross-coupling error matrix.
5155 * @param m23 element 2,3 of cross-coupling error matrix.
5156 * @param m33 element 3,3 of cross-coupling error matrix.
5157 * @param g11 element 1,1 of g-dependent cross bias matrix.
5158 * @param g21 element 2,1 of g-dependent cross bias matrix.
5159 * @param g31 element 3,1 of g-dependent cross bias matrix.
5160 * @param g12 element 1,2 of g-dependent cross bias matrix.
5161 * @param g22 element 2,2 of g-dependent cross bias matrix.
5162 * @param g32 element 3,2 of g-dependent cross bias matrix.
5163 * @param g13 element 1,3 of g-dependent cross bias matrix.
5164 * @param g23 element 2,3 of g-dependent cross bias matrix.
5165 * @param g33 element 3,3 of g-dependent cross bias matrix.
5166 * @return error between estimated and measured gravity versor.
5167 * @throws EvaluationException if there are numerical instabilities.
5168 */
5169 private double evaluate(
5170 final int i, final double bx, final double by, final double bz,
5171 final double m11, final double m21, final double m31,
5172 final double m12, final double m22, final double m32,
5173 final double m13, final double m23, final double m33,
5174 final double g11, final double g21, final double g31,
5175 final double g12, final double g22, final double g32,
5176 final double g13, final double g23, final double g33) throws EvaluationException {
5177
5178 try {
5179 final var measuredSequence = sequences.get(i);
5180 final var fixedSequence = fixedSequences.get(i);
5181
5182 // generate new sequence using current parameters to fix angular rate measurements
5183 if (measuredSpecificForce == null) {
5184 measuredSpecificForce = new Matrix(BodyKinematics.COMPONENTS, 1);
5185 }
5186 if (trueSpecificForce == null) {
5187 trueSpecificForce = new Matrix(BodyKinematics.COMPONENTS, 1);
5188 }
5189
5190 if (measuredAngularRate == null) {
5191 measuredAngularRate = new Matrix(BodyKinematics.COMPONENTS, 1);
5192 }
5193 if (trueAngularRate == null) {
5194 trueAngularRate = new Matrix(BodyKinematics.COMPONENTS, 1);
5195 }
5196
5197 if (m == null) {
5198 m = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
5199 }
5200 if (invM == null) {
5201 invM = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
5202 }
5203 if (b == null) {
5204 b = new Matrix(BodyKinematics.COMPONENTS, 1);
5205 }
5206 if (g == null) {
5207 g = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
5208 }
5209 if (tmp == null) {
5210 tmp = new Matrix(BodyKinematics.COMPONENTS, 1);
5211 }
5212
5213 m.setElementAt(0, 0, m11);
5214 m.setElementAt(1, 0, m21);
5215 m.setElementAt(2, 0, m31);
5216
5217 m.setElementAt(0, 1, m12);
5218 m.setElementAt(1, 1, m22);
5219 m.setElementAt(2, 1, m32);
5220
5221 m.setElementAt(0, 2, m13);
5222 m.setElementAt(1, 2, m23);
5223 m.setElementAt(2, 2, m33);
5224
5225 Utils.inverse(m, invM);
5226
5227 b.setElementAtIndex(0, bx);
5228 b.setElementAtIndex(1, by);
5229 b.setElementAtIndex(2, bz);
5230
5231 g.setElementAt(0, 0, g11);
5232 g.setElementAt(1, 0, g21);
5233 g.setElementAt(2, 0, g31);
5234
5235 g.setElementAt(0, 1, g12);
5236 g.setElementAt(1, 1, g22);
5237 g.setElementAt(2, 1, g32);
5238
5239 g.setElementAt(0, 2, g13);
5240 g.setElementAt(1, 2, g23);
5241 g.setElementAt(2, 2, g33);
5242
5243 // fix kinematics
5244 final var numItems = measuredSequence.getItemsCount();
5245 final var measuredItems = measuredSequence.getSortedItems();
5246 final var fixedItems = fixedSequence.getSortedItems();
5247 for (var j = 0; j < numItems; j++) {
5248 final var measuredItem = measuredItems.get(j);
5249 final var fixedItem = fixedItems.get(j);
5250
5251 fixKinematics(measuredItem.getKinematics(), fixedItem.getKinematics());
5252 }
5253
5254
5255 // integrate fixed sequence to obtain attitude change
5256 QuaternionIntegrator.integrateGyroSequence(fixedSequence, QuaternionStepIntegratorType.RUNGE_KUTTA, q);
5257
5258 startPoint.setInhomogeneousCoordinates(point[0], point[1], point[2]);
5259 q.inverse();
5260 q.rotate(startPoint, endPoint);
5261
5262 expectedEndPoint.setInhomogeneousCoordinates(point[3], point[4], point[5]);
5263
5264 return expectedEndPoint.distanceTo(endPoint);
5265
5266 } catch (final AlgebraException | RotationException e) {
5267 throw new EvaluationException(e);
5268 }
5269 }
5270
5271 /**
5272 * Fixes provided kinematics with provided accelerometer parameters and
5273 * current gyroscope parameters.
5274 *
5275 * @param kinematics kinematics to be fixed with current values.
5276 * @param result kinematics where result will be stored.
5277 * @throws AlgebraException if for some reason kinematics
5278 */
5279 private void fixKinematics(final BodyKinematics kinematics, final BodyKinematics result) throws AlgebraException {
5280
5281 // Ωmeas = bg + (I + Mg) * Ωtrue + Gg * ftrue
5282 // Ωmeas = M*(Ωtrue + b + G * ftrue)
5283
5284 // M = I + Mg
5285 // bg = M*b --> b = M^-1*bg
5286 // Gg = M*G --> G = M^-1*Gg
5287
5288 // Ωtrue = M^-1 * Ωmeas - b - G*ftrue
5289
5290 // fix specific force
5291 measuredSpecificForce.setElementAtIndex(0, kinematics.getFx());
5292 measuredSpecificForce.setElementAtIndex(1, kinematics.getFy());
5293 measuredSpecificForce.setElementAtIndex(2, kinematics.getFz());
5294
5295 accelerationFixer.fix(measuredSpecificForce, trueSpecificForce);
5296
5297 // fix angular rate
5298 measuredAngularRate.setElementAtIndex(0, kinematics.getAngularRateX());
5299 measuredAngularRate.setElementAtIndex(1, kinematics.getAngularRateY());
5300 measuredAngularRate.setElementAtIndex(2, kinematics.getAngularRateZ());
5301
5302 g.multiply(trueSpecificForce, tmp);
5303 invM.multiply(measuredAngularRate, trueAngularRate);
5304 trueAngularRate.subtract(b);
5305 trueAngularRate.subtract(tmp);
5306
5307 result.setSpecificForceCoordinates(
5308 trueSpecificForce.getElementAtIndex(0),
5309 trueSpecificForce.getElementAtIndex(1),
5310 trueSpecificForce.getElementAtIndex(2));
5311
5312 result.setAngularRateCoordinates(
5313 trueAngularRate.getElementAtIndex(0),
5314 trueAngularRate.getElementAtIndex(1),
5315 trueAngularRate.getElementAtIndex(2));
5316 }
5317
5318 /**
5319 * Resets any previous estimations.
5320 */
5321 private void reset() {
5322 estimatedBiases = null;
5323 estimatedMg = null;
5324 estimatedGg = null;
5325 estimatedCovariance = null;
5326 estimatedChiSq = 0.0;
5327 estimatedChiSqDegreesOfFreedom = 0;
5328 estimatedReducedChiSq = 0.0;
5329 estimatedMse = 0.0;
5330 estimatedP = 1.0;
5331 estimatedQ = 0.0;
5332 i = 0;
5333 fixedSequences = null;
5334 measuredSpecificForce = null;
5335 trueSpecificForce = null;
5336 measuredAngularRate = null;
5337 trueAngularRate = null;
5338 m = null;
5339 invM = null;
5340 b = null;
5341 g = null;
5342 tmp = null;
5343 point = null;
5344 }
5345 }