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