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.Matrix;
20  import com.irurueta.algebra.Utils;
21  import com.irurueta.algebra.WrongSizeException;
22  import com.irurueta.navigation.LockedException;
23  import com.irurueta.navigation.NotReadyException;
24  import com.irurueta.navigation.frames.CoordinateTransformation;
25  import com.irurueta.navigation.frames.ECEFPosition;
26  import com.irurueta.navigation.frames.ECEFVelocity;
27  import com.irurueta.navigation.frames.FrameType;
28  import com.irurueta.navigation.frames.InvalidSourceAndDestinationFrameTypeException;
29  import com.irurueta.navigation.frames.NEDFrame;
30  import com.irurueta.navigation.frames.NEDPosition;
31  import com.irurueta.navigation.frames.NEDVelocity;
32  import com.irurueta.navigation.frames.converters.ECEFtoNEDPositionVelocityConverter;
33  import com.irurueta.navigation.frames.converters.NEDtoECEFFrameConverter;
34  import com.irurueta.navigation.frames.converters.NEDtoECEFPositionVelocityConverter;
35  import com.irurueta.navigation.geodesic.Constants;
36  import com.irurueta.navigation.inertial.BodyKinematics;
37  import com.irurueta.navigation.inertial.INSLooselyCoupledKalmanInitializerConfig;
38  import com.irurueta.navigation.inertial.INSTightlyCoupledKalmanInitializerConfig;
39  import com.irurueta.navigation.inertial.calibration.AccelerationFixer;
40  import com.irurueta.navigation.inertial.calibration.AngularSpeedTriad;
41  import com.irurueta.navigation.inertial.calibration.CalibrationException;
42  import com.irurueta.navigation.inertial.calibration.GyroscopeBiasUncertaintySource;
43  import com.irurueta.navigation.inertial.calibration.GyroscopeCalibrationSource;
44  import com.irurueta.navigation.inertial.calibration.StandardDeviationBodyKinematics;
45  import com.irurueta.navigation.inertial.estimators.ECEFKinematicsEstimator;
46  import com.irurueta.numerical.EvaluationException;
47  import com.irurueta.numerical.GradientEstimator;
48  import com.irurueta.numerical.fitting.FittingException;
49  import com.irurueta.numerical.fitting.LevenbergMarquardtMultiDimensionFitter;
50  import com.irurueta.numerical.fitting.LevenbergMarquardtMultiDimensionFunctionEvaluator;
51  import com.irurueta.statistics.MaxIterationsExceededException;
52  import com.irurueta.units.Acceleration;
53  import com.irurueta.units.AccelerationConverter;
54  import com.irurueta.units.AccelerationUnit;
55  import com.irurueta.units.AngularSpeed;
56  import com.irurueta.units.AngularSpeedConverter;
57  import com.irurueta.units.AngularSpeedUnit;
58  import com.irurueta.units.Time;
59  import com.irurueta.units.TimeConverter;
60  import com.irurueta.units.TimeUnit;
61  
62  import java.util.Collection;
63  
64  /**
65   * Estimates gyroscope biases, cross couplings and scaling factors
66   * along with G-dependent cross biases introduced on the gyroscope by the
67   * specific forces sensed by the accelerometer.
68   * <p>
69   * This calibrator assumes that the IMU is placed flat on a turntable spinning
70   * at constant speed, but absolute orientation or position of IMU is unknown.
71   * Turntable must rotate fast enough so that Earth rotation effects can be
72   * neglected, bus slow enough so that gyroscope readings can be properly made.
73   * <p>
74   * To use this calibrator at least 10 measurements are needed when common
75   * z-axis is assumed and G-dependent cross biases are ignored, otherwise
76   * at least 13 measurements are required when common z-axis is not assumed.
77   * If G-dependent cross biases are being estimated, then at least 19
78   * measurements are needed when common z-axis is assumed, otherwise at
79   * least 22 measurements are required when common z-axis is not assumed.
80   * <p>
81   * Measured gyroscope angular rates is assumed to follow the model shown below:
82   * <pre>
83   *     Ωmeas = bg + (I + Mg) * Ωtrue + Gg * ftrue + w
84   * </pre>
85   * Where:
86   * - Ωmeas is the measured gyroscope angular rates. This is a 3x1 vector.
87   * - bg is the gyroscope bias. Ideally, on a perfect gyroscope, this should be a
88   * 3x1 zero vector.
89   * - I is the 3x3 identity matrix.
90   * - Mg is the 3x3 matrix containing cross-couplings and scaling factors. Ideally, on
91   * a perfect gyroscope, this should be a 3x3 zero matrix.
92   * - Ωtrue is ground-truth gyroscope angular rates.
93   * - Gg is the G-dependent cross biases introduced by the specific forces sensed
94   * by the accelerometer. Ideally, on a perfect gyroscope, this should be a 3x3
95   * zero matrix.
96   * - ftrue is ground-truth specific force. This is a 3x1 vector.
97   * - w is measurement noise. This is a 3x1 vector.
98   */
99  public class TurntableGyroscopeCalibrator implements GyroscopeNonLinearCalibrator, UnknownBiasGyroscopeCalibrator,
100         GyroscopeCalibrationSource, GyroscopeBiasUncertaintySource,
101         UnorderedStandardDeviationBodyKinematicsGyroscopeCalibrator, AccelerometerDependentGyroscopeCalibrator {
102 
103     /**
104      * Indicates whether by default a common z-axis is assumed for both the accelerometer
105      * and gyroscope.
106      */
107     public static final boolean DEFAULT_USE_COMMON_Z_AXIS = true;
108 
109     /**
110      * Indicates that by default G-dependent cross biases introduced
111      * by the accelerometer on the gyroscope are estimated.
112      */
113     public static final boolean DEFAULT_ESTIMATE_G_DEPENDENT_CROSS_BIASES = true;
114 
115     /**
116      * Number of unknowns when common z-axis is assumed for both the accelerometer
117      * and gyroscope when G-dependent cross biases are being estimated.
118      */
119     public static final int COMMON_Z_AXIS_UNKNOWNS_AND_CROSS_BIASES = 18;
120 
121     /**
122      * Number of unknowns for the general case when G-dependent cross
123      * biases are being estimated.
124      */
125     public static final int GENERAL_UNKNOWNS_AND_CROSS_BIASES = 21;
126 
127     /**
128      * Number of unknowns when common z-axis is assumed for both
129      * the accelerometer and gyroscope when G-dependent cross biases
130      * are not being estimated.
131      */
132     public static final int COMMON_Z_AXIS_UNKNOWNS = 9;
133 
134     /**
135      * Number of unknowns for the general case when G-dependent cross
136      * biases are not being estimated.
137      */
138     public static final int GENERAL_UNKNOWNS = 12;
139 
140     /**
141      * Required minimum number of measurements when common z-axis is assumed
142      * and G-dependent cross biases are being estimated.
143      */
144     public static final int MINIMUM_MEASUREMENTS_COMMON_Z_AXIS_AND_CROSS_BIASES =
145             COMMON_Z_AXIS_UNKNOWNS_AND_CROSS_BIASES + 1;
146 
147     /**
148      * Required minimum number of measurements for the general case and
149      * G-dependent cross biases are being estimated.
150      */
151     public static final int MINIMUM_MEASUREMENTS_GENERAL_AND_CROSS_BIASES = GENERAL_UNKNOWNS_AND_CROSS_BIASES + 1;
152 
153     /**
154      * Required minimum number of measurements when common z-axis is assumed
155      * and G-dependent cross biases are being ignored.
156      */
157     public static final int MINIMUM_MEASUREMENTS_COMMON_Z_AXIS = COMMON_Z_AXIS_UNKNOWNS + 1;
158 
159     /**
160      * Required minimum number of measurements for the general case and
161      * G-dependent cross biases are being ignored.
162      */
163     public static final int MINIMUM_MEASUREMENTS_GENERAL = GENERAL_UNKNOWNS + 1;
164 
165     /**
166      * Default turntable rotation rate.
167      */
168     public static final double DEFAULT_TURNTABLE_ROTATION_RATE = Constants.EARTH_ROTATION_RATE;
169 
170     /**
171      * Default time interval between measurements expressed in seconds (s).
172      * This is a typical value when we have 50 samples per second.
173      */
174     public static final double DEFAULT_TIME_INTERVAL = 0.02;
175 
176     /**
177      * Levenberg-Marquardt fitter to find a non-linear solution.
178      */
179     private final LevenbergMarquardtMultiDimensionFitter fitter = new LevenbergMarquardtMultiDimensionFitter();
180 
181     /**
182      * Known x-coordinate of accelerometer bias to be used to fix measured
183      * specific force and find cross biases introduced by the accelerometer.
184      * This is expressed in meters per squared second (m/s^2).
185      */
186     private double accelerometerBiasX;
187 
188     /**
189      * Known y-coordinate of accelerometer bias to be used to fix measured
190      * specific force and find cross biases introduced by the accelerometer.
191      * This is expressed in meters per squared second (m/s^2).
192      */
193     private double accelerometerBiasY;
194 
195     /**
196      * Known z-coordinate of accelerometer bias to be used to fix measured
197      * specific force and find cross biases introduced by the accelerometer.
198      * This is expressed in meters per squared second (m/s^2).
199      */
200     private double accelerometerBiasZ;
201 
202     /**
203      * Known accelerometer x scaling factor to be used to fix measured
204      * specific force and find cross biases introduced by the accelerometer.
205      */
206     private double accelerometerSx;
207 
208     /**
209      * Known accelerometer y scaling factor to be used to fix measured
210      * specific force and find cross biases introduced by the accelerometer.
211      */
212     private double accelerometerSy;
213 
214     /**
215      * Known accelerometer z scaling factor to be used to fix measured
216      * specific force and find cross biases introduced by the accelerometer.
217      */
218     private double accelerometerSz;
219 
220     /**
221      * Known accelerometer x-y cross coupling error to be used to fix measured
222      * specific force and find cross biases introduced by the accelerometer.
223      */
224     private double accelerometerMxy;
225 
226     /**
227      * Know accelerometer x-z cross coupling error to be used to fix measured
228      * specific force and find cross biases introduced by the accelerometer.
229      */
230     private double accelerometerMxz;
231 
232     /**
233      * Known accelerometer y-x cross coupling error to be used to fix measured
234      * specific force and find cross biases introduced by the accelerometer.
235      */
236     private double accelerometerMyx;
237 
238     /**
239      * Known accelerometer y-z cross coupling error to be used to fix measured
240      * specific force and find cross biases introduced by the accelerometer.
241      */
242     private double accelerometerMyz;
243 
244     /**
245      * Known accelerometer z-x cross coupling error to be used to fix measured
246      * specific force and find cross biases introduced by the accelerometer.
247      */
248     private double accelerometerMzx;
249 
250     /**
251      * Known accelerometer z-y cross coupling error to be used to fix measured
252      * specific force and find cross biases introduced by the accelerometer.
253      */
254     private double accelerometerMzy;
255 
256     /**
257      * Initial x-coordinate of gyroscope bias to be used to find a solution.
258      * This is expressed in radians per second (rad/s).
259      */
260     private double initialBiasX;
261 
262     /**
263      * Initial y-coordinate of gyroscope bias to be used to find a solution.
264      * This is expressed in radians per second (rad/s).
265      */
266     private double initialBiasY;
267 
268     /**
269      * Initial z-coordinate of gyroscope bias to be used to find a solution.
270      * This is expressed in radians per second (rad/s).
271      */
272     private double initialBiasZ;
273 
274     /**
275      * Initial gyroscope x scaling factor.
276      */
277     private double initialSx;
278 
279     /**
280      * Initial gyroscope y scaling factor.
281      */
282     private double initialSy;
283 
284     /**
285      * Initial gyroscope z scaling factor.
286      */
287     private double initialSz;
288 
289     /**
290      * Initial gyroscope x-y cross coupling error.
291      */
292     private double initialMxy;
293 
294     /**
295      * Initial gyroscope x-z cross coupling error.
296      */
297     private double initialMxz;
298 
299     /**
300      * Initial gyroscope y-x cross coupling error.
301      */
302     private double initialMyx;
303 
304     /**
305      * Initial gyroscope y-z cross coupling error.
306      */
307     private double initialMyz;
308 
309     /**
310      * Initial gyroscope z-x cross coupling error.
311      */
312     private double initialMzx;
313 
314     /**
315      * Initial gyroscope z-y cross coupling error.
316      */
317     private double initialMzy;
318 
319     /**
320      * Initial G-dependent cross biases introduced on the gyroscope by the
321      * specific forces sensed by the accelerometer.
322      */
323     private Matrix initialGg;
324 
325     /**
326      * Constant rotation rate at which the turntable is spinning.
327      * This is expressed in radians per second (rad/s).
328      */
329     private double turntableRotationRate = DEFAULT_TURNTABLE_ROTATION_RATE;
330 
331     /**
332      * Time interval between measurements being captured expressed in
333      * second (s).
334      */
335     private double timeInterval = DEFAULT_TIME_INTERVAL;
336 
337     /**
338      * Contains a collection of body kinematics measurements taken at
339      * a given position with different unknown orientations and containing
340      * the standard deviations of accelerometer and gyroscope measurements.
341      */
342     private Collection<StandardDeviationBodyKinematics> measurements;
343 
344     /**
345      * Position where body kinematics measures have been taken.
346      */
347     private ECEFPosition position;
348 
349     /**
350      * This flag indicates whether z-axis is assumed to be common for accelerometer
351      * and gyroscope.
352      * When enabled, this eliminates 3 variables from Mg matrix.
353      */
354     private boolean commonAxisUsed = DEFAULT_USE_COMMON_Z_AXIS;
355 
356     /**
357      * This flag indicates whether G-dependent cross biases are being
358      * estimated or not.
359      * When enabled, this adds 9 variables from Gg matrix.
360      */
361     private boolean estimateGDependentCrossBiases = DEFAULT_ESTIMATE_G_DEPENDENT_CROSS_BIASES;
362 
363     /**
364      * Listener to handle events raised by this calibrator.
365      */
366     private TurntableGyroscopeCalibratorListener listener;
367 
368     /**
369      * Estimated angular rate biases for each IMU axis expressed in radians per
370      * second (rad/s).
371      */
372     private double[] estimatedBiases;
373 
374     /**
375      * Estimated gyroscope scale factors and cross coupling errors.
376      * This is the product of matrix Tg containing cross coupling errors and Kg
377      * containing scaling factors.
378      * So that:
379      * <pre>
380      *     Mg = [sx    mxy  mxz] = Tg*Kg
381      *          [myx   sy   myz]
382      *          [mzx   mzy  sz ]
383      * </pre>
384      * Where:
385      * <pre>
386      *     Kg = [sx 0   0 ]
387      *          [0  sy  0 ]
388      *          [0  0   sz]
389      * </pre>
390      * and
391      * <pre>
392      *     Tg = [1          -alphaXy    alphaXz ]
393      *          [alphaYx    1           -alphaYz]
394      *          [-alphaZx   alphaZy     1       ]
395      * </pre>
396      * Hence:
397      * <pre>
398      *     Mg = [sx    mxy  mxz] = Tg*Kg =  [sx             -sy * alphaXy   sz * alphaXz ]
399      *          [myx   sy   myz]            [sx * alphaYx   sy              -sz * alphaYz]
400      *          [mzx   mzy  sz ]            [-sx * alphaZx  sy * alphaZy    sz           ]
401      * </pre>
402      * This instance allows any 3x3 matrix however, typically alphaYx, alphaZx and alphaZy
403      * are considered to be zero if the gyroscope z-axis is assumed to be the same
404      * as the body z-axis. When this is assumed, myx = mzx = mzy = 0 and the Mg matrix
405      * becomes upper diagonal:
406      * <pre>
407      *     Mg = [sx    mxy  mxz]
408      *          [0     sy   myz]
409      *          [0     0    sz ]
410      * </pre>
411      * Values of this matrix are unit-less.
412      */
413     private Matrix estimatedMg;
414 
415     /**
416      * Estimated G-dependent cross biases introduced on the gyroscope by the
417      * specific forces sensed by the accelerometer.
418      * This instance allows any 3x3 matrix.
419      */
420     private Matrix estimatedGg;
421 
422     /**
423      * Estimated covariance matrix for estimated parameters.
424      */
425     private Matrix estimatedCovariance;
426 
427     /**
428      * Estimated chi square value.
429      */
430     private double estimatedChiSq;
431 
432     /**
433      * Estimated degrees of freedom of chi square value. Degrees of freedom is equal to the number of sampled data
434      * minus the number of estimated parameters.
435      */
436     private int estimatedChiSqDegreesOfFreedom;
437 
438     /**
439      * Estimated reduced chi square value. This is equal to estimated chi square value divided by its degrees of
440      * freedom. Ideally this value should be close to 1.0.
441      */
442     private double estimatedReducedChiSq;
443 
444     /**
445      * Estimated mean square error respect to provided measurements.
446      */
447     private double estimatedMse;
448 
449     /**
450      * Estimated probability of finding a smaller chi square value expressed as a value between 0.0 and 1.0. The smaller
451      * the found chi square value is, the better the fit of the estimated parameters to the actual parameter. Thus, the
452      * smaller the chance of finding a smaller chi square value, then the better the estimated fit is.
453      */
454     private double estimatedP;
455 
456     /**
457      * Estimated measure of quality of estimated fit as a value between 0.0 and 1.0. The larger the quality value is,
458      * the better the fit that has been estimated.
459      */
460     private double estimatedQ;
461 
462     /**
463      * Indicates whether calibrator is running.
464      */
465     private boolean running;
466 
467     /**
468      * Internally holds x-coordinate of measured angular rate during calibration.
469      */
470     private double measAngularRateX;
471 
472     /**
473      * Internally holds y-coordinate of measured angular rate during calibration.
474      */
475     private double measAngularRateY;
476 
477     /**
478      * Internally holds z-coordinate of measured angular rate during calibration.
479      */
480     private double measAngularRateZ;
481 
482     /**
483      * Internally holds x-coordinate of measured specific force during calibration.
484      */
485     private double fmeasX;
486 
487     /**
488      * Internally holds y-coordinate of measured specific force during calibration.
489      */
490     private double fmeasY;
491 
492     /**
493      * Internally holds z-coordinate of measured specific force during calibration.
494      */
495     private double fmeasZ;
496 
497     /**
498      * Internally holds measured angular rate during calibration expressed as
499      * a column matrix.
500      */
501     private Matrix measAngularRate;
502 
503     /**
504      * Internally holds measured specific force during calibration expressed as
505      * a column matrix.
506      */
507     private Matrix fmeas;
508 
509     /**
510      * Internally holds cross-coupling errors during calibration.
511      */
512     private Matrix m;
513 
514     /**
515      * Internally holds inverse of cross-coupling errors during calibration.
516      */
517     private Matrix invM;
518 
519     /**
520      * Internally holds biases during calibration.
521      */
522     private Matrix b;
523 
524     /**
525      * Internally hold g-dependent cross biases during calibration.
526      */
527     private Matrix g;
528 
529     /**
530      * Internally holds computed true angular rate during calibration.
531      */
532     private Matrix trueAngularRate;
533 
534     /**
535      * Internally holds computed true specific force during calibration.
536      */
537     private Matrix ftrue;
538 
539     /**
540      * Internally holds accelerometer bias during calibration.
541      */
542     private Matrix ba;
543 
544     /**
545      * Internally holds accelerometer scaling and cross coupling errors
546      * during calibration.
547      */
548     private Matrix ma;
549 
550     /**
551      * Internally holds angular rate bias due to g-dependent cross biases
552      */
553     private Matrix tmp;
554 
555     /**
556      * Acceleration fixer.
557      */
558     private final AccelerationFixer accelerationFixer = new AccelerationFixer();
559 
560     /**
561      * Constructor.
562      */
563     public TurntableGyroscopeCalibrator() {
564         try {
565             initialGg = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
566         } catch (final WrongSizeException ignore) {
567             // never happens
568         }
569     }
570 
571     /**
572      * Constructor.
573      *
574      * @param position              position where body kinematics measures
575      *                              have been taken.
576      * @param turntableRotationRate constant rotation rate at which the
577      *                              turntable is spinning. Must be
578      *                              expressed in radians per second (rad/s).
579      * @param timeInterval          time interval between measurements being
580      *                              captured expressed in seconds (s).
581      * @param measurements          collection of body kinematics
582      *                              measurements with standard deviations
583      *                              taken at the same position with zero
584      *                              velocity and unknown different
585      *                              orientations.
586      * @param initialBias           initial gyroscope bias to be used to
587      *                              find a solution. This must be 3x1 and
588      *                              is expressed in radians per second
589      *                              (rad/s).
590      * @param initialMg             initial gyroscope scale factors and
591      *                              cross coupling errors matrix. Must
592      *                              be 3x3.
593      * @param initialGg             initial gyroscope G-dependent cross
594      *                              biases introduced on the gyroscope by
595      *                              the specific forces sensed by the
596      *                              accelerometer. Must be 3x3.
597      * @throws IllegalArgumentException if any of the provided values does
598      *                                  not have proper size or if either
599      *                                  turntable rotation rate or
600      *                                  time interval is zero or negative.
601      */
602     public TurntableGyroscopeCalibrator(
603             final ECEFPosition position, final double turntableRotationRate, final double timeInterval,
604             final Collection<StandardDeviationBodyKinematics> measurements, final Matrix initialBias,
605             final Matrix initialMg, final Matrix initialGg) {
606         this();
607         this.position = position;
608         this.measurements = measurements;
609         try {
610             setTurntableRotationRate(turntableRotationRate);
611             setTimeInterval(timeInterval);
612             setInitialBias(initialBias);
613             setInitialMg(initialMg);
614             setInitialGg(initialGg);
615         } catch (final LockedException ignore) {
616             // never happens
617         }
618     }
619 
620     /**
621      * Constructor.
622      *
623      * @param position              position where body kinematics measures
624      *                              have been taken.
625      * @param turntableRotationRate constant rotation rate at which the
626      *                              turntable is spinning. Must be
627      *                              expressed in radians per second (rad/s).
628      * @param timeInterval          time interval between measurements being
629      *                              captured expressed in seconds (s).
630      * @param measurements          collection of body kinematics
631      *                              measurements with standard deviations
632      *                              taken at the same position with zero
633      *                              velocity and unknown different
634      *                              orientations.
635      * @param initialBias           initial gyroscope bias to be used to
636      *                              find a solution. This must be 3x1 and
637      *                              is expressed in radians per second
638      *                              (rad/s).
639      * @param initialMg             initial gyroscope scale factors and
640      *                              cross coupling errors matrix. Must
641      *                              be 3x3.
642      * @param initialGg             initial gyroscope G-dependent cross
643      *                              biases introduced on the gyroscope by
644      *                              the specific forces sensed by the
645      *                              accelerometer. Must be 3x3.
646      * @param listener              listener to handle events raised by this
647      *                              calibrator.
648      * @throws IllegalArgumentException if any of the provided values does
649      *                                  not have proper size or if either
650      *                                  turntable rotation rate or
651      *                                  time interval is zero or negative.
652      */
653     public TurntableGyroscopeCalibrator(
654             final ECEFPosition position, final double turntableRotationRate, final double timeInterval,
655             final Collection<StandardDeviationBodyKinematics> measurements, final Matrix initialBias,
656             final Matrix initialMg, final Matrix initialGg, final TurntableGyroscopeCalibratorListener listener) {
657         this(position, turntableRotationRate, timeInterval, measurements, initialBias, initialMg, initialGg);
658         this.listener = listener;
659     }
660 
661     /**
662      * Constructor.
663      *
664      * @param position              position where body kinematics measures
665      *                              have been taken.
666      * @param turntableRotationRate constant rotation rate at which the
667      *                              turntable is spinning. Must be
668      *                              expressed in radians per second (rad/s).
669      * @param timeInterval          time interval between measurements being
670      *                              captured expressed in seconds (s).
671      * @param measurements          collection of body kinematics
672      *                              measurements with standard deviations
673      *                              taken at the same position with zero
674      *                              velocity and unknown different
675      *                              orientations.
676      * @param initialBias           initial gyroscope bias to be used to
677      *                              find a solution. This must have
678      *                              length 3 and is expressed in radians
679      *                              per second (rad/s).
680      * @param initialMg             initial gyroscope scale factors and
681      *                              cross coupling errors matrix. Must
682      *                              be 3x3.
683      * @param initialGg             initial gyroscope G-dependent cross
684      *                              biases introduced on the gyroscope by
685      *                              the specific forces sensed by the
686      *                              accelerometer. Must be 3x3.
687      * @throws IllegalArgumentException if any of the provided values does
688      *                                  not have proper size or if either
689      *                                  turntable rotation rate or
690      *                                  time interval is zero or negative.
691      */
692     public TurntableGyroscopeCalibrator(
693             final ECEFPosition position, final double turntableRotationRate, final double timeInterval,
694             final Collection<StandardDeviationBodyKinematics> measurements, final double[] initialBias,
695             final Matrix initialMg, final Matrix initialGg) {
696         this();
697         this.position = position;
698         this.measurements = measurements;
699         try {
700             setTurntableRotationRate(turntableRotationRate);
701             setTimeInterval(timeInterval);
702             setInitialBias(initialBias);
703             setInitialMg(initialMg);
704             setInitialGg(initialGg);
705         } catch (final LockedException ignore) {
706             // never happens
707         }
708     }
709 
710     /**
711      * Constructor.
712      *
713      * @param position              position where body kinematics measures
714      *                              have been taken.
715      * @param turntableRotationRate constant rotation rate at which the
716      *                              turntable is spinning. Must be
717      *                              expressed in radians per second (rad/s).
718      * @param timeInterval          time interval between measurements being
719      *                              captured expressed in seconds (s).
720      * @param measurements          collection of body kinematics
721      *                              measurements with standard deviations
722      *                              taken at the same position with zero
723      *                              velocity and unknown different
724      *                              orientations.
725      * @param initialBias           initial gyroscope bias to be used to
726      *                              find a solution. This must have length
727      *                              3 and is expressed in radians
728      *                              per second (rad/s).
729      * @param initialMg             initial gyroscope scale factors and
730      *                              cross coupling errors matrix. Must
731      *                              be 3x3.
732      * @param initialGg             initial gyroscope G-dependent cross
733      *                              biases introduced on the gyroscope by
734      *                              the specific forces sensed by the
735      *                              accelerometer. Must be 3x3.
736      * @param listener              listener to handle events raised by
737      *                              this calibrator.
738      * @throws IllegalArgumentException if any of the provided values does
739      *                                  not have proper size or if either
740      *                                  turntable rotation rate or
741      *                                  time interval is zero or negative.
742      */
743     public TurntableGyroscopeCalibrator(
744             final ECEFPosition position, final double turntableRotationRate, final double timeInterval,
745             final Collection<StandardDeviationBodyKinematics> measurements, final double[] initialBias,
746             final Matrix initialMg, final Matrix initialGg, final TurntableGyroscopeCalibratorListener listener) {
747         this(position, turntableRotationRate, timeInterval, measurements, initialBias, initialMg, initialGg);
748         this.listener = listener;
749     }
750 
751     /**
752      * Constructor.
753      *
754      * @param position              position where body kinematics measures
755      *                              have been taken.
756      * @param turntableRotationRate constant rotation rate at which the
757      *                              turntable is spinning. Must be
758      *                              expressed in radians per second (rad/s).
759      * @param timeInterval          time interval between measurements being
760      *                              captured expressed in seconds (s).
761      * @param measurements          collection of body kinematics
762      *                              measurements with standard deviations
763      *                              taken at the same position with zero
764      *                              velocity and unknown different
765      *                              orientations.
766      * @param initialBias           initial gyroscope bias to be used to
767      *                              find a solution. This must have length
768      *                              3 and is expressed in radians per
769      *                              second (rad/s).
770      * @param initialMg             initial gyroscope scale factors and
771      *                              cross coupling errors matrix. Must
772      *                              be 3x3.
773      * @param initialGg             initial gyroscope G-dependent cross
774      *                              biases introduced on the gyroscope by
775      *                              the specific forces sensed by the
776      *                              accelerometer. Must be 3x3.
777      * @param accelerometerBias     known accelerometer bias. This must
778      *                              have length 3 and is expressed in
779      *                              meters per squared second
780      *                              (m/s^2).
781      * @param accelerometerMa       known accelerometer scale factors and
782      *                              cross coupling matrix. Must be 3x3.
783      * @throws IllegalArgumentException if any of the provided values does
784      *                                  not have proper size or if either
785      *                                  turntable rotation rate or
786      *                                  time interval is zero or negative.
787      */
788     public TurntableGyroscopeCalibrator(
789             final ECEFPosition position, final double turntableRotationRate, final double timeInterval,
790             final Collection<StandardDeviationBodyKinematics> measurements, final double[] initialBias,
791             final Matrix initialMg, final Matrix initialGg, final double[] accelerometerBias,
792             final Matrix accelerometerMa) {
793         this(position, turntableRotationRate, timeInterval, measurements, initialBias, initialMg, initialGg);
794         try {
795             setAccelerometerBias(accelerometerBias);
796             setAccelerometerMa(accelerometerMa);
797         } catch (final LockedException ignore) {
798             // never happens
799         }
800     }
801 
802     /**
803      * Constructor.
804      *
805      * @param position              position where body kinematics measures
806      *                              have been taken.
807      * @param turntableRotationRate constant rotation rate at which the
808      *                              turntable is spinning. Must be
809      *                              expressed in radians per second (rad/s).
810      * @param timeInterval          time interval between measurements being
811      *                              captured expressed in seconds (s).
812      * @param measurements          collection of body kinematics
813      *                              measurements with standard deviations
814      *                              taken at the same position with zero
815      *                              velocity and unknown different
816      *                              orientations.
817      * @param initialBias           initial gyroscope bias to be used to
818      *                              find a solution. This must have length
819      *                              3 and is expressed in radians per
820      *                              second (rad/s).
821      * @param initialMg             initial gyroscope scale factors and
822      *                              cross coupling errors matrix. Must
823      *                              be 3x3.
824      * @param initialGg             initial gyroscope G-dependent cross
825      *                              biases introduced on the gyroscope by
826      *                              the specific forces sensed by the
827      *                              accelerometer. Must be 3x3.
828      * @param accelerometerBias     known accelerometer bias. This must
829      *                              have length 3 and is expressed in
830      *                              meters per squared second (m/s^2).
831      * @param accelerometerMa       known accelerometer scale factors and
832      *                              cross coupling matrix. Must be 3x3.
833      * @param listener              listener to handle events raised by
834      *                              this calibrator.
835      * @throws IllegalArgumentException if any of the provided values does
836      *                                  not have proper size or if either
837      *                                  turntable rotation rate or
838      *                                  time interval is zero or negative.
839      */
840     public TurntableGyroscopeCalibrator(
841             final ECEFPosition position, final double turntableRotationRate, final double timeInterval,
842             final Collection<StandardDeviationBodyKinematics> measurements, final double[] initialBias,
843             final Matrix initialMg, final Matrix initialGg, final double[] accelerometerBias,
844             final Matrix accelerometerMa, final TurntableGyroscopeCalibratorListener listener) {
845         this(position, turntableRotationRate, timeInterval, measurements, initialBias, initialMg, initialGg,
846                 accelerometerBias, accelerometerMa);
847         this.listener = listener;
848     }
849 
850     /**
851      * Constructor.
852      *
853      * @param position              position where body kinematics measures
854      *                              have been taken.
855      * @param turntableRotationRate constant rotation rate at which the
856      *                              turntable is spinning. Must be
857      *                              expressed in radians per second (rad/s).
858      * @param timeInterval          time interval between measurements being
859      *                              captured expressed in seconds (s).
860      * @param measurements          collection of body kinematics
861      *                              measurements with standard deviations
862      *                              taken at the same position with zero
863      *                              velocity and unknown different
864      *                              orientations.
865      * @param initialBias           initial gyroscope bias to be used to
866      *                              find a solution. This must be 3x1 and
867      *                              is expressed in radians per second
868      *                              (rad/s).
869      * @param initialMg             initial gyroscope scale factors and
870      *                              cross coupling errors matrix. Must
871      *                              be 3x3.
872      * @param initialGg             initial gyroscope G-dependent cross
873      *                              biases introduced on the gyroscope by
874      *                              the specific forces sensed by the
875      *                              accelerometer. Must be 3x3.
876      * @param accelerometerBias     known accelerometer bias. This must be 3x1
877      *                              and is expressed in meters per squared
878      *                              second (m/s^2).
879      * @param accelerometerMa       known accelerometer scale factors and
880      *                              cross coupling matrix. Must be 3x3.
881      * @throws IllegalArgumentException if any of the provided values does
882      *                                  not have proper size or if either
883      *                                  turntable rotation rate or
884      *                                  time interval is zero or negative.
885      */
886     public TurntableGyroscopeCalibrator(
887             final ECEFPosition position, final double turntableRotationRate, final double timeInterval,
888             final Collection<StandardDeviationBodyKinematics> measurements, final Matrix initialBias,
889             final Matrix initialMg, final Matrix initialGg, final Matrix accelerometerBias,
890             final Matrix accelerometerMa) {
891         this(position, turntableRotationRate, timeInterval, measurements, initialBias, initialMg, initialGg);
892         try {
893             setAccelerometerBias(accelerometerBias);
894             setAccelerometerMa(accelerometerMa);
895         } catch (final LockedException ignore) {
896             // never happens
897         }
898     }
899 
900     /**
901      * Constructor.
902      *
903      * @param position              position where body kinematics measures
904      *                              have been taken.
905      * @param turntableRotationRate constant rotation rate at which the
906      *                              turntable is spinning. Must be
907      *                              expressed in radians per second (rad/s).
908      * @param timeInterval          time interval between measurements being
909      *                              captured expressed in seconds (s).
910      * @param measurements          collection of body kinematics
911      *                              measurements with standard deviations
912      *                              taken at the same position with zero
913      *                              velocity and unknown different
914      *                              orientations.
915      * @param initialBias           initial gyroscope bias to be used to
916      *                              find a solution. This must be 3x1 and
917      *                              is expressed in radians per second
918      *                              (rad/s).
919      * @param initialMg             initial gyroscope scale factors and
920      *                              cross coupling errors matrix. Must
921      *                              be 3x3.
922      * @param initialGg             initial gyroscope G-dependent cross
923      *                              biases introduced on the gyroscope by
924      *                              the specific forces sensed by the
925      *                              accelerometer. Must be 3x3.
926      * @param accelerometerBias     known accelerometer bias. This must
927      *                              have length 3 and is expressed in
928      *                              meters per squared second (m/s^2).
929      * @param accelerometerMa       known accelerometer scale factors and
930      *                              cross coupling matrix. Must be 3x3.
931      * @param listener              listener to handle events raised by
932      *                              this calibrator.
933      * @throws IllegalArgumentException if any of the provided values does
934      *                                  not have proper size or if either
935      *                                  turntable rotation rate or
936      *                                  time interval is zero or negative.
937      */
938     public TurntableGyroscopeCalibrator(
939             final ECEFPosition position, final double turntableRotationRate, final double timeInterval,
940             final Collection<StandardDeviationBodyKinematics> measurements, final Matrix initialBias,
941             final Matrix initialMg, final Matrix initialGg, final Matrix accelerometerBias,
942             final Matrix accelerometerMa, final TurntableGyroscopeCalibratorListener listener) {
943         this(position, turntableRotationRate, timeInterval, measurements, initialBias, initialMg, initialGg,
944                 accelerometerBias, accelerometerMa);
945         this.listener = listener;
946     }
947 
948     /**
949      * Constructor.
950      *
951      * @param position                      position where body kinematics
952      *                                      measures have been taken.
953      * @param turntableRotationRate         constant rotation rate at which
954      *                                      the turntable is spinning. Must
955      *                                      be expressed in radians per
956      *                                      second (rad/s).
957      * @param timeInterval                  time interval between measurements
958      *                                      being captured expressed in
959      *                                      seconds (s).
960      * @param measurements                  collection of body kinematics
961      *                                      measurements with standard
962      *                                      deviations taken at the same
963      *                                      position with zero velocity
964      *                                      and unknown different
965      *                                      orientations.
966      * @param commonAxisUsed                indicates whether z-axis is
967      *                                      assumed to be common for
968      *                                      accelerometer and gyroscope.
969      * @param estimateGDependentCrossBiases true if G-dependent cross biases
970      *                                      will be estimated, false
971      *                                      otherwise.
972      * @param initialBias                   initial gyroscope bias to be
973      *                                      used to find a solution. This
974      *                                      must be 3x1 and is expressed in
975      *                                      radians per second (rad/s).
976      * @param initialMg                     initial gyroscope scale factors
977      *                                      and cross coupling errors matrix.
978      *                                      Must be 3x3.
979      * @param initialGg                     initial gyroscope G-dependent
980      *                                      cross biases introduced on the
981      *                                      gyroscope by the specific
982      *                                      forces sensed by the
983      *                                      accelerometer. Must be 3x3.
984      * @throws IllegalArgumentException if any of the provided values does
985      *                                  not have proper size or if either
986      *                                  turntable rotation rate or
987      *                                  time interval is zero or negative.
988      */
989     public TurntableGyroscopeCalibrator(
990             final ECEFPosition position, final double turntableRotationRate, final double timeInterval,
991             final Collection<StandardDeviationBodyKinematics> measurements, final boolean commonAxisUsed,
992             final boolean estimateGDependentCrossBiases, final Matrix initialBias, final Matrix initialMg,
993             final Matrix initialGg) {
994         this(position, turntableRotationRate, timeInterval, measurements, initialBias, initialMg, initialGg);
995         this.commonAxisUsed = commonAxisUsed;
996         this.estimateGDependentCrossBiases = estimateGDependentCrossBiases;
997     }
998 
999     /**
1000      * Constructor.
1001      *
1002      * @param position                      position where body kinematics
1003      *                                      measures have been taken.
1004      * @param turntableRotationRate         constant rotation rate at which
1005      *                                      the turntable is spinning. Must
1006      *                                      be expressed in radians per
1007      *                                      second (rad/s).
1008      * @param timeInterval                  time interval between measurements
1009      *                                      being captured expressed in
1010      *                                      seconds (s).
1011      * @param measurements                  collection of body kinematics
1012      *                                      measurements with standard
1013      *                                      deviations taken at the same
1014      *                                      position with zero velocity and
1015      *                                      unknown different orientations.
1016      * @param commonAxisUsed                indicates whether z-axis is
1017      *                                      assumed to be common for
1018      *                                      accelerometer and gyroscope.
1019      * @param estimateGDependentCrossBiases true if G-dependent cross
1020      *                                      biases will be estimated, false
1021      *                                      otherwise.
1022      * @param initialBias                   initial gyroscope bias to be
1023      *                                      used to find a solution. This
1024      *                                      must be 3x1 and is expressed in
1025      *                                      radians per second (rad/s).
1026      * @param initialMg                     initial gyroscope scale factors
1027      *                                      and cross coupling errors
1028      *                                      matrix. Must be 3x3.
1029      * @param initialGg                     initial gyroscope G-dependent
1030      *                                      cross biases introduced on the
1031      *                                      gyroscope by the specific
1032      *                                      forces sensed by the
1033      *                                      accelerometer. Must be 3x3.
1034      * @param listener                      listener to handle events
1035      *                                      raised by this calibrator.
1036      * @throws IllegalArgumentException if any of the provided values does
1037      *                                  not have proper size or if either
1038      *                                  turntable rotation rate or
1039      *                                  time interval is zero or negative.
1040      */
1041     public TurntableGyroscopeCalibrator(
1042             final ECEFPosition position, final double turntableRotationRate, final double timeInterval,
1043             final Collection<StandardDeviationBodyKinematics> measurements, final boolean commonAxisUsed,
1044             final boolean estimateGDependentCrossBiases, final Matrix initialBias, final Matrix initialMg,
1045             final Matrix initialGg, final TurntableGyroscopeCalibratorListener listener) {
1046         this(position, turntableRotationRate, timeInterval, measurements, commonAxisUsed,
1047                 estimateGDependentCrossBiases, initialBias, initialMg, initialGg);
1048         this.listener = listener;
1049     }
1050 
1051     /**
1052      * Constructor.
1053      *
1054      * @param position                      position where body kinematics
1055      *                                      measures have been taken.
1056      * @param turntableRotationRate         constant rotation rate at which
1057      *                                      the turntable is spinning. Must
1058      *                                      be expressed in radians per
1059      *                                      second (rad/s).
1060      * @param timeInterval                  time interval between measurements
1061      *                                      being captured expressed in
1062      *                                      seconds (s).
1063      * @param measurements                  collection of body kinematics
1064      *                                      measurements with standard
1065      *                                      deviations taken at the same
1066      *                                      position with zero velocity
1067      *                                      and unknown different
1068      *                                      orientations.
1069      * @param commonAxisUsed                indicates whether z-axis is
1070      *                                      assumed to be common for
1071      *                                      accelerometer and gyroscope.
1072      * @param estimateGDependentCrossBiases true if G-dependent cross biases
1073      *                                      will be estimated, false
1074      *                                      otherwise.
1075      * @param initialBias                   initial gyroscope bias to be
1076      *                                      used to find a solution. This
1077      *                                      must have length 3 and is
1078      *                                      expressed in radians per second
1079      *                                      (rad/s).
1080      * @param initialMg                     initial gyroscope scale factors
1081      *                                      and cross coupling errors matrix.
1082      *                                      Must be 3x3.
1083      * @param initialGg                     initial gyroscope G-dependent
1084      *                                      cross biases introduced on the
1085      *                                      gyroscope by the specific forces
1086      *                                      sensed by the accelerometer.
1087      *                                      Must be 3x3.
1088      * @throws IllegalArgumentException if any of the provided values does
1089      *                                  not have proper size or if either
1090      *                                  turntable rotation rate or
1091      *                                  time interval is zero or negative.
1092      */
1093     public TurntableGyroscopeCalibrator(
1094             final ECEFPosition position, final double turntableRotationRate, final double timeInterval,
1095             final Collection<StandardDeviationBodyKinematics> measurements, final boolean commonAxisUsed,
1096             final boolean estimateGDependentCrossBiases, final double[] initialBias, final Matrix initialMg,
1097             final Matrix initialGg) {
1098         this(position, turntableRotationRate, timeInterval, measurements, initialBias, initialMg, initialGg);
1099         this.commonAxisUsed = commonAxisUsed;
1100         this.estimateGDependentCrossBiases = estimateGDependentCrossBiases;
1101     }
1102 
1103     /**
1104      * Constructor.
1105      *
1106      * @param position                      position where body kinematics
1107      *                                      measures have been taken.
1108      * @param turntableRotationRate         constant rotation rate at which
1109      *                                      the turntable is spinning. Must
1110      *                                      be expressed in radians per
1111      *                                      second (rad/s).
1112      * @param timeInterval                  time interval between measurements
1113      *                                      being captured expressed in
1114      *                                      seconds (s).
1115      * @param measurements                  collection of body kinematics
1116      *                                      measurements with standard
1117      *                                      deviations taken at the same
1118      *                                      position with zero velocity
1119      *                                      and unknown different
1120      *                                      orientations.
1121      * @param commonAxisUsed                indicates whether z-axis is
1122      *                                      assumed to be common for
1123      *                                      accelerometer and gyroscope.
1124      * @param estimateGDependentCrossBiases true if G-dependent cross biases
1125      *                                      will be estimated, false
1126      *                                      otherwise.
1127      * @param initialBias                   initial gyroscope bias to be
1128      *                                      used to find a solution. This
1129      *                                      must have length 3 and is
1130      *                                      expressed in radians per second
1131      *                                      (rad/s).
1132      * @param initialMg                     initial gyroscope scale factors
1133      *                                      and cross coupling errors
1134      *                                      matrix. Must be 3x3.
1135      * @param initialGg                     initial gyroscope G-dependent
1136      *                                      cross biases introduced on the
1137      *                                      gyroscope by the specific forces
1138      *                                      sensed by the accelerometer.
1139      *                                      Must be 3x3.
1140      * @param listener                      listener to handle events raised
1141      *                                      by this calibrator.
1142      * @throws IllegalArgumentException if any of the provided values does
1143      *                                  not have proper size or if either
1144      *                                  turntable rotation rate or
1145      *                                  time interval is zero or negative.
1146      */
1147     public TurntableGyroscopeCalibrator(
1148             final ECEFPosition position, final double turntableRotationRate, final double timeInterval,
1149             final Collection<StandardDeviationBodyKinematics> measurements, final boolean commonAxisUsed,
1150             final boolean estimateGDependentCrossBiases, final double[] initialBias, final Matrix initialMg,
1151             final Matrix initialGg, final TurntableGyroscopeCalibratorListener listener) {
1152         this(position, turntableRotationRate, timeInterval, measurements, commonAxisUsed, estimateGDependentCrossBiases,
1153                 initialBias, initialMg, initialGg);
1154         this.listener = listener;
1155     }
1156 
1157     /**
1158      * Constructor.
1159      *
1160      * @param position                      position where body kinematics
1161      *                                      measures have been taken.
1162      * @param turntableRotationRate         constant rotation rate at which
1163      *                                      the turntable is spinning. Must
1164      *                                      be expressed in radians per
1165      *                                      second (rad/s).
1166      * @param timeInterval                  time interval between measurements
1167      *                                      being captured expressed in
1168      *                                      seconds (s).
1169      * @param measurements                  collection of body kinematics
1170      *                                      measurements with standard
1171      *                                      deviations taken at the same
1172      *                                      position with zero velocity
1173      *                                      and unknown different
1174      *                                      orientations.
1175      * @param commonAxisUsed                indicates whether z-axis is
1176      *                                      assumed to be common for
1177      *                                      accelerometer and gyroscope.
1178      * @param estimateGDependentCrossBiases true if G-dependent cross
1179      *                                      biases will be estimated,
1180      *                                      false otherwise.
1181      * @param initialBias                   initial gyroscope bias to be
1182      *                                      used to find a solution. This
1183      *                                      must have length 3 and is
1184      *                                      expressed in radians per second
1185      *                                      (rad/s).
1186      * @param initialMg                     initial gyroscope scale factors
1187      *                                      and cross coupling errors
1188      *                                      matrix. Must be 3x3.
1189      * @param initialGg                     initial gyroscope G-dependent
1190      *                                      cross biases introduced on the
1191      *                                      gyroscope by the specific forces
1192      *                                      sensed by the accelerometer.
1193      *                                      Must be 3x3.
1194      * @param accelerometerBias             known accelerometer bias. This
1195      *                                      must have length 3 and is
1196      *                                      expressed in meters per squared
1197      *                                      second (m/s^2).
1198      * @param accelerometerMa               known accelerometer scale factors
1199      *                                      and cross coupling matrix. Must
1200      *                                      be 3x3.
1201      * @throws IllegalArgumentException if any of the provided values does
1202      *                                  not have proper size or if either
1203      *                                  turntable rotation rate or
1204      *                                  time interval is zero or negative.
1205      */
1206     public TurntableGyroscopeCalibrator(
1207             final ECEFPosition position, final double turntableRotationRate, final double timeInterval,
1208             final Collection<StandardDeviationBodyKinematics> measurements, final boolean commonAxisUsed,
1209             final boolean estimateGDependentCrossBiases, final double[] initialBias, final Matrix initialMg,
1210             final Matrix initialGg, final double[] accelerometerBias, final Matrix accelerometerMa) {
1211         this(position, turntableRotationRate, timeInterval, measurements, initialBias, initialMg, initialGg,
1212                 accelerometerBias, accelerometerMa);
1213         this.commonAxisUsed = commonAxisUsed;
1214         this.estimateGDependentCrossBiases = estimateGDependentCrossBiases;
1215     }
1216 
1217     /**
1218      * Constructor.
1219      *
1220      * @param position                      position where body kinematics
1221      *                                      measures have been taken.
1222      * @param turntableRotationRate         constant rotation rate at which
1223      *                                      the turntable is spinning. Must
1224      *                                      be expressed in radians per
1225      *                                      second (rad/s).
1226      * @param timeInterval                  time interval between measurements
1227      *                                      being captured expressed in
1228      *                                      seconds (s).
1229      * @param measurements                  collection of body kinematics
1230      *                                      measurements with standard
1231      *                                      deviations taken at the same
1232      *                                      position with zero velocity
1233      *                                      and unknown different
1234      *                                      orientations.
1235      * @param commonAxisUsed                indicates whether z-axis is
1236      *                                      assumed to be common for
1237      *                                      accelerometer and gyroscope.
1238      * @param estimateGDependentCrossBiases true if G-dependent cross biases
1239      *                                      will be estimated, false
1240      *                                      otherwise.
1241      * @param initialBias                   initial gyroscope bias to be used
1242      *                                      to find a solution. This must
1243      *                                      have length 3 and is expressed
1244      *                                      in radians per second (rad/s).
1245      * @param initialMg                     initial gyroscope scale factors
1246      *                                      and cross coupling errors matrix.
1247      *                                      Must be 3x3.
1248      * @param initialGg                     initial gyroscope G-dependent
1249      *                                      cross biases introduced on the
1250      *                                      gyroscope by the specific forces
1251      *                                      sensed by the accelerometer. Must
1252      *                                      be 3x3.
1253      * @param accelerometerBias             known accelerometer bias. This
1254      *                                      must have length 3 and is
1255      *                                      expressed in meters per squared
1256      *                                      second (m/s^2).
1257      * @param accelerometerMa               known accelerometer scale factors
1258      *                                      and cross coupling matrix. Must
1259      *                                      be 3x3.
1260      * @param listener                      listener to handle events raised
1261      *                                      by this calibrator.
1262      * @throws IllegalArgumentException if any of the provided values does
1263      *                                  not have proper size or if either
1264      *                                  turntable rotation rate or
1265      *                                  time interval is zero or negative.
1266      */
1267     public TurntableGyroscopeCalibrator(
1268             final ECEFPosition position, final double turntableRotationRate, final double timeInterval,
1269             final Collection<StandardDeviationBodyKinematics> measurements, final boolean commonAxisUsed,
1270             final boolean estimateGDependentCrossBiases, final double[] initialBias, final Matrix initialMg,
1271             final Matrix initialGg, final double[] accelerometerBias, final Matrix accelerometerMa,
1272             final TurntableGyroscopeCalibratorListener listener) {
1273         this(position, turntableRotationRate, timeInterval, measurements, commonAxisUsed, estimateGDependentCrossBiases,
1274                 initialBias, initialMg, initialGg, accelerometerBias, accelerometerMa);
1275         this.listener = listener;
1276     }
1277 
1278     /**
1279      * Constructor.
1280      *
1281      * @param position                      position where body kinematics
1282      *                                      measures have been taken.
1283      * @param turntableRotationRate         constant rotation rate at which
1284      *                                      the turntable is spinning. Must
1285      *                                      be expressed in radians per
1286      *                                      second (rad/s).
1287      * @param timeInterval                  time interval between measurements
1288      *                                      being captured expressed in
1289      *                                      seconds (s).
1290      * @param measurements                  collection of body kinematics
1291      *                                      measurements with standard
1292      *                                      deviations taken at the same
1293      *                                      position with zero velocity and
1294      *                                      unknown different orientations.
1295      * @param commonAxisUsed                indicates whether z-axis is
1296      *                                      assumed to be common for
1297      *                                      accelerometer and gyroscope.
1298      * @param estimateGDependentCrossBiases true if G-dependent cross biases
1299      *                                      will be estimated, false
1300      *                                      otherwise.
1301      * @param initialBias                   initial gyroscope bias to be
1302      *                                      used to find a solution. This
1303      *                                      must be 3x1 and is expressed in
1304      *                                      radians per second (rad/s).
1305      * @param initialMg                     initial gyroscope scale factors
1306      *                                      and cross coupling errors matrix.
1307      *                                      Must be 3x3.
1308      * @param initialGg                     initial gyroscope G-dependent
1309      *                                      cross biases introduced on the
1310      *                                      gyroscope by the specific forces
1311      *                                      sensed by the accelerometer. Must
1312      *                                      be 3x3.
1313      * @param accelerometerBias             known accelerometer bias. This
1314      *                                      must have length 3 and is
1315      *                                      expressed in meters per squared
1316      *                                      second (m/s^2).
1317      * @param accelerometerMa               known accelerometer scale factors
1318      *                                      and cross coupling matrix. Must
1319      *                                      be 3x3.
1320      * @throws IllegalArgumentException if any of the provided values does
1321      *                                  not have proper size or if either
1322      *                                  turntable rotation rate or
1323      *                                  time interval is zero or negative.
1324      */
1325     public TurntableGyroscopeCalibrator(
1326             final ECEFPosition position, final double turntableRotationRate, final double timeInterval,
1327             final Collection<StandardDeviationBodyKinematics> measurements, final boolean commonAxisUsed,
1328             final boolean estimateGDependentCrossBiases, final Matrix initialBias, final Matrix initialMg,
1329             final Matrix initialGg, final Matrix accelerometerBias, final Matrix accelerometerMa) {
1330         this(position, turntableRotationRate, timeInterval, measurements, initialBias, initialMg, initialGg,
1331                 accelerometerBias, accelerometerMa);
1332         this.commonAxisUsed = commonAxisUsed;
1333         this.estimateGDependentCrossBiases = estimateGDependentCrossBiases;
1334     }
1335 
1336     /**
1337      * Constructor.
1338      *
1339      * @param position                      position where body kinematics
1340      *                                      measures have been taken.
1341      * @param turntableRotationRate         constant rotation rate at which
1342      *                                      the turntable is spinning. Must
1343      *                                      be expressed in radians per
1344      *                                      second (rad/s).
1345      * @param timeInterval                  time interval between measurements
1346      *                                      being captured expressed in
1347      *                                      seconds (s).
1348      * @param measurements                  collection of body kinematics
1349      *                                      measurements with standard
1350      *                                      deviations taken at the same
1351      *                                      position with zero velocity and
1352      *                                      unknown different orientations.
1353      * @param commonAxisUsed                indicates whether z-axis is
1354      *                                      assumed to be common for
1355      *                                      accelerometer and gyroscope.
1356      * @param estimateGDependentCrossBiases true if G-dependent cross biases
1357      *                                      will be estimated, false
1358      *                                      otherwise.
1359      * @param initialBias                   initial gyroscope bias to be used
1360      *                                      to find a solution. This must be
1361      *                                      3x1 and is expressed in radians
1362      *                                      per second (rad/s).
1363      * @param initialMg                     initial gyroscope scale factors
1364      *                                      and cross coupling errors matrix.
1365      *                                      Must be 3x3.
1366      * @param initialGg                     initial gyroscope G-dependent
1367      *                                      cross biases introduced on the
1368      *                                      gyroscope by the specific forces
1369      *                                      sensed by the accelerometer. Must
1370      *                                      be 3x3.
1371      * @param accelerometerBias             known accelerometer bias. This
1372      *                                      must have length 3 and is
1373      *                                      expressed in meters per squared
1374      *                                      second (m/s^2).
1375      * @param accelerometerMa               known accelerometer scale factors
1376      *                                      and cross coupling matrix. Must
1377      *                                      be 3x3.
1378      * @param listener                      listener to handle events raised
1379      *                                      by this calibrator.
1380      * @throws IllegalArgumentException if any of the provided values does
1381      *                                  not have proper size or if either
1382      *                                  turntable rotation rate or
1383      *                                  time interval is zero or negative.
1384      */
1385     public TurntableGyroscopeCalibrator(
1386             final ECEFPosition position, final double turntableRotationRate, final double timeInterval,
1387             final Collection<StandardDeviationBodyKinematics> measurements, final boolean commonAxisUsed,
1388             final boolean estimateGDependentCrossBiases, final Matrix initialBias, final Matrix initialMg,
1389             final Matrix initialGg, final Matrix accelerometerBias, final Matrix accelerometerMa,
1390             final TurntableGyroscopeCalibratorListener listener) {
1391         this(position, turntableRotationRate, timeInterval, measurements, commonAxisUsed,
1392                 estimateGDependentCrossBiases, initialBias, initialMg, initialGg, accelerometerBias, accelerometerMa);
1393         this.listener = listener;
1394     }
1395 
1396     /**
1397      * Constructor.
1398      *
1399      * @param position              position where body kinematics measures
1400      *                              have been taken.
1401      * @param turntableRotationRate constant rotation rate at which the
1402      *                              turntable is spinning. Must be
1403      *                              expressed in radians per second (rad/s).
1404      * @param timeInterval          time interval between measurements being
1405      *                              captured expressed in seconds (s).
1406      * @param measurements          collection of body kinematics
1407      *                              measurements with standard deviations
1408      *                              taken at the same position with zero
1409      *                              velocity and unknown different
1410      *                              orientations.
1411      * @param initialBias           initial gyroscope bias to be used to
1412      *                              find a solution. This must be 3x1 and
1413      *                              is expressed in radians per second
1414      *                              (rad/s).
1415      * @param initialMg             initial gyroscope scale factors and
1416      *                              cross coupling errors matrix. Must
1417      *                              be 3x3.
1418      * @param initialGg             initial gyroscope G-dependent cross
1419      *                              biases introduced on the gyroscope by
1420      *                              the specific forces sensed by the
1421      *                              accelerometer. Must be 3x3.
1422      * @throws IllegalArgumentException if any of the provided values does
1423      *                                  not have proper size or if either
1424      *                                  turntable rotation rate or
1425      *                                  time interval is zero or negative.
1426      */
1427     public TurntableGyroscopeCalibrator(
1428             final NEDPosition position, final double turntableRotationRate, final double timeInterval,
1429             final Collection<StandardDeviationBodyKinematics> measurements, final Matrix initialBias,
1430             final Matrix initialMg, final Matrix initialGg) {
1431         this(convertPosition(position), turntableRotationRate, timeInterval, measurements, initialBias, initialMg,
1432                 initialGg);
1433     }
1434 
1435     /**
1436      * Constructor.
1437      *
1438      * @param position              position where body kinematics measures
1439      *                              have been taken.
1440      * @param turntableRotationRate constant rotation rate at which the
1441      *                              turntable is spinning. Must be
1442      *                              expressed in radians per second (rad/s).
1443      * @param timeInterval          time interval between measurements being
1444      *                              captured expressed in seconds (s).
1445      * @param measurements          collection of body kinematics
1446      *                              measurements with standard deviations
1447      *                              taken at the same position with zero
1448      *                              velocity and unknown different
1449      *                              orientations.
1450      * @param initialBias           initial gyroscope bias to be used to
1451      *                              find a solution. This must be 3x1 and
1452      *                              is expressed in radians per second
1453      *                              (rad/s).
1454      * @param initialMg             initial gyroscope scale factors and
1455      *                              cross coupling errors matrix. Must
1456      *                              be 3x3.
1457      * @param initialGg             initial gyroscope G-dependent cross
1458      *                              biases introduced on the gyroscope by
1459      *                              the specific forces sensed by the
1460      *                              accelerometer. Must be 3x3.
1461      * @param listener              listener to handle events raised by this
1462      *                              calibrator.
1463      * @throws IllegalArgumentException if any of the provided values does
1464      *                                  not have proper size or if either
1465      *                                  turntable rotation rate or
1466      *                                  time interval is zero or negative.
1467      */
1468     public TurntableGyroscopeCalibrator(
1469             final NEDPosition position, final double turntableRotationRate, final double timeInterval,
1470             final Collection<StandardDeviationBodyKinematics> measurements, final Matrix initialBias,
1471             final Matrix initialMg, final Matrix initialGg, final TurntableGyroscopeCalibratorListener listener) {
1472         this(convertPosition(position), turntableRotationRate, timeInterval, measurements, initialBias, initialMg,
1473                 initialGg, listener);
1474     }
1475 
1476     /**
1477      * Constructor.
1478      *
1479      * @param position              position where body kinematics measures
1480      *                              have been taken.
1481      * @param turntableRotationRate constant rotation rate at which the
1482      *                              turntable is spinning. Must be
1483      *                              expressed in radians per second (rad/s).
1484      * @param timeInterval          time interval between measurements being
1485      *                              captured expressed in seconds (s).
1486      * @param measurements          collection of body kinematics
1487      *                              measurements with standard deviations
1488      *                              taken at the same position with zero
1489      *                              velocity and unknown different
1490      *                              orientations.
1491      * @param initialBias           initial gyroscope bias to be used to
1492      *                              find a solution. This must have
1493      *                              length 3 and is expressed in radians
1494      *                              per second (rad/s).
1495      * @param initialMg             initial gyroscope scale factors and
1496      *                              cross coupling errors matrix. Must
1497      *                              be 3x3.
1498      * @param initialGg             initial gyroscope G-dependent cross
1499      *                              biases introduced on the gyroscope by
1500      *                              the specific forces sensed by the
1501      *                              accelerometer. Must be 3x3.
1502      * @throws IllegalArgumentException if any of the provided values does
1503      *                                  not have proper size or if either
1504      *                                  turntable rotation rate or
1505      *                                  time interval is zero or negative.
1506      */
1507     public TurntableGyroscopeCalibrator(
1508             final NEDPosition position, final double turntableRotationRate, final double timeInterval,
1509             final Collection<StandardDeviationBodyKinematics> measurements, final double[] initialBias,
1510             final Matrix initialMg, final Matrix initialGg) {
1511         this(convertPosition(position), turntableRotationRate, timeInterval, measurements, initialBias, initialMg,
1512                 initialGg);
1513     }
1514 
1515     /**
1516      * Constructor.
1517      *
1518      * @param position              position where body kinematics measures
1519      *                              have been taken.
1520      * @param turntableRotationRate constant rotation rate at which the
1521      *                              turntable is spinning. Must be
1522      *                              expressed in radians per second (rad/s).
1523      * @param timeInterval          time interval between measurements being
1524      *                              captured expressed in seconds (s).
1525      * @param measurements          collection of body kinematics
1526      *                              measurements with standard deviations
1527      *                              taken at the same position with zero
1528      *                              velocity and unknown different
1529      *                              orientations.
1530      * @param initialBias           initial gyroscope bias to be used to
1531      *                              find a solution. This must have length
1532      *                              3 and is expressed in radians
1533      *                              per second (rad/s).
1534      * @param initialMg             initial gyroscope scale factors and
1535      *                              cross coupling errors matrix. Must
1536      *                              be 3x3.
1537      * @param initialGg             initial gyroscope G-dependent cross
1538      *                              biases introduced on the gyroscope by
1539      *                              the specific forces sensed by the
1540      *                              accelerometer. Must be 3x3.
1541      * @param listener              listener to handle events raised by
1542      *                              this calibrator.
1543      * @throws IllegalArgumentException if any of the provided values does
1544      *                                  not have proper size or if either
1545      *                                  turntable rotation rate or
1546      *                                  time interval is zero or negative.
1547      */
1548     public TurntableGyroscopeCalibrator(
1549             final NEDPosition position, final double turntableRotationRate, final double timeInterval,
1550             final Collection<StandardDeviationBodyKinematics> measurements, final double[] initialBias,
1551             final Matrix initialMg, final Matrix initialGg, final TurntableGyroscopeCalibratorListener listener) {
1552         this(convertPosition(position), turntableRotationRate, timeInterval, measurements, initialBias, initialMg,
1553                 initialGg, listener);
1554     }
1555 
1556     /**
1557      * Constructor.
1558      *
1559      * @param position              position where body kinematics measures
1560      *                              have been taken.
1561      * @param turntableRotationRate constant rotation rate at which the
1562      *                              turntable is spinning. Must be
1563      *                              expressed in radians per second (rad/s).
1564      * @param timeInterval          time interval between measurements being
1565      *                              captured expressed in seconds (s).
1566      * @param measurements          collection of body kinematics
1567      *                              measurements with standard deviations
1568      *                              taken at the same position with zero
1569      *                              velocity and unknown different
1570      *                              orientations.
1571      * @param initialBias           initial gyroscope bias to be used to
1572      *                              find a solution. This must have length
1573      *                              3 and is expressed in radians per
1574      *                              second (rad/s).
1575      * @param initialMg             initial gyroscope scale factors and
1576      *                              cross coupling errors matrix. Must
1577      *                              be 3x3.
1578      * @param initialGg             initial gyroscope G-dependent cross
1579      *                              biases introduced on the gyroscope by
1580      *                              the specific forces sensed by the
1581      *                              accelerometer. Must be 3x3.
1582      * @param accelerometerBias     known accelerometer bias. This must
1583      *                              have length 3 and is expressed in
1584      *                              meters per squared second
1585      *                              (m/s^2).
1586      * @param accelerometerMa       known accelerometer scale factors and
1587      *                              cross coupling matrix. Must be 3x3.
1588      * @throws IllegalArgumentException if any of the provided values does
1589      *                                  not have proper size or if either
1590      *                                  turntable rotation rate or
1591      *                                  time interval is zero or negative.
1592      */
1593     public TurntableGyroscopeCalibrator(
1594             final NEDPosition position, final double turntableRotationRate, final double timeInterval,
1595             final Collection<StandardDeviationBodyKinematics> measurements, final double[] initialBias,
1596             final Matrix initialMg, final Matrix initialGg, final double[] accelerometerBias,
1597             final Matrix accelerometerMa) {
1598         this(convertPosition(position), turntableRotationRate, timeInterval, measurements, initialBias, initialMg,
1599                 initialGg, accelerometerBias, accelerometerMa);
1600     }
1601 
1602     /**
1603      * Constructor.
1604      *
1605      * @param position              position where body kinematics measures
1606      *                              have been taken.
1607      * @param turntableRotationRate constant rotation rate at which the
1608      *                              turntable is spinning. Must be
1609      *                              expressed in radians per second (rad/s).
1610      * @param timeInterval          time interval between measurements being
1611      *                              captured expressed in seconds (s).
1612      * @param measurements          collection of body kinematics
1613      *                              measurements with standard deviations
1614      *                              taken at the same position with zero
1615      *                              velocity and unknown different
1616      *                              orientations.
1617      * @param initialBias           initial gyroscope bias to be used to
1618      *                              find a solution. This must have length
1619      *                              3 and is expressed in radians per
1620      *                              second (rad/s).
1621      * @param initialMg             initial gyroscope scale factors and
1622      *                              cross coupling errors matrix. Must
1623      *                              be 3x3.
1624      * @param initialGg             initial gyroscope G-dependent cross
1625      *                              biases introduced on the gyroscope by
1626      *                              the specific forces sensed by the
1627      *                              accelerometer. Must be 3x3.
1628      * @param accelerometerBias     known accelerometer bias. This must
1629      *                              have length 3 and is expressed in
1630      *                              meters per squared second (m/s^2).
1631      * @param accelerometerMa       known accelerometer scale factors and
1632      *                              cross coupling matrix. Must be 3x3.
1633      * @param listener              listener to handle events raised by
1634      *                              this calibrator.
1635      * @throws IllegalArgumentException if any of the provided values does
1636      *                                  not have proper size or if either
1637      *                                  turntable rotation rate or
1638      *                                  time interval is zero or negative.
1639      */
1640     public TurntableGyroscopeCalibrator(
1641             final NEDPosition position, final double turntableRotationRate, final double timeInterval,
1642             final Collection<StandardDeviationBodyKinematics> measurements, final double[] initialBias,
1643             final Matrix initialMg, final Matrix initialGg, final double[] accelerometerBias,
1644             final Matrix accelerometerMa, final TurntableGyroscopeCalibratorListener listener) {
1645         this(convertPosition(position), turntableRotationRate, timeInterval, measurements, initialBias, initialMg,
1646                 initialGg, accelerometerBias, accelerometerMa, listener);
1647     }
1648 
1649     /**
1650      * Constructor.
1651      *
1652      * @param position              position where body kinematics measures
1653      *                              have been taken.
1654      * @param turntableRotationRate constant rotation rate at which the
1655      *                              turntable is spinning. Must be
1656      *                              expressed in radians per second (rad/s).
1657      * @param timeInterval          time interval between measurements being
1658      *                              captured expressed in seconds (s).
1659      * @param measurements          collection of body kinematics
1660      *                              measurements with standard deviations
1661      *                              taken at the same position with zero
1662      *                              velocity and unknown different
1663      *                              orientations.
1664      * @param initialBias           initial gyroscope bias to be used to
1665      *                              find a solution. This must be 3x1 and
1666      *                              is expressed in radians per second
1667      *                              (rad/s).
1668      * @param initialMg             initial gyroscope scale factors and
1669      *                              cross coupling errors matrix. Must
1670      *                              be 3x3.
1671      * @param initialGg             initial gyroscope G-dependent cross
1672      *                              biases introduced on the gyroscope by
1673      *                              the specific forces sensed by the
1674      *                              accelerometer. Must be 3x3.
1675      * @param accelerometerBias     known accelerometer bias. This must
1676      *                              have length 3 and is expressed in
1677      *                              meters per squared second
1678      *                              (m/s^2).
1679      * @param accelerometerMa       known accelerometer scale factors and
1680      *                              cross coupling matrix. Must be 3x3.
1681      * @throws IllegalArgumentException if any of the provided values does
1682      *                                  not have proper size or if either
1683      *                                  turntable rotation rate or
1684      *                                  time interval is zero or negative.
1685      */
1686     public TurntableGyroscopeCalibrator(
1687             final NEDPosition position, final double turntableRotationRate, final double timeInterval,
1688             final Collection<StandardDeviationBodyKinematics> measurements, final Matrix initialBias,
1689             final Matrix initialMg, final Matrix initialGg, final Matrix accelerometerBias,
1690             final Matrix accelerometerMa) {
1691         this(convertPosition(position), turntableRotationRate, timeInterval, measurements, initialBias, initialMg,
1692                 initialGg, accelerometerBias, accelerometerMa);
1693     }
1694 
1695     /**
1696      * Constructor.
1697      *
1698      * @param position              position where body kinematics measures
1699      *                              have been taken.
1700      * @param turntableRotationRate constant rotation rate at which the
1701      *                              turntable is spinning. Must be
1702      *                              expressed in radians per second (rad/s).
1703      * @param timeInterval          time interval between measurements being
1704      *                              captured expressed in seconds (s).
1705      * @param measurements          collection of body kinematics
1706      *                              measurements with standard deviations
1707      *                              taken at the same position with zero
1708      *                              velocity and unknown different
1709      *                              orientations.
1710      * @param initialBias           initial gyroscope bias to be used to
1711      *                              find a solution. This must be 3x1 and
1712      *                              is expressed in radians per second
1713      *                              (rad/s).
1714      * @param initialMg             initial gyroscope scale factors and
1715      *                              cross coupling errors matrix. Must
1716      *                              be 3x3.
1717      * @param initialGg             initial gyroscope G-dependent cross
1718      *                              biases introduced on the gyroscope by
1719      *                              the specific forces sensed by the
1720      *                              accelerometer. Must be 3x3.
1721      * @param accelerometerBias     known accelerometer bias. This must
1722      *                              have length 3 and is expressed in
1723      *                              meters per squared second (m/s^2).
1724      * @param accelerometerMa       known accelerometer scale factors and
1725      *                              cross coupling matrix. Must be 3x3.
1726      * @param listener              listener to handle events raised by
1727      *                              this calibrator.
1728      * @throws IllegalArgumentException if any of the provided values does
1729      *                                  not have proper size or if either
1730      *                                  turntable rotation rate or
1731      *                                  time interval is zero or negative.
1732      */
1733     public TurntableGyroscopeCalibrator(
1734             final NEDPosition position, final double turntableRotationRate, final double timeInterval,
1735             final Collection<StandardDeviationBodyKinematics> measurements, final Matrix initialBias,
1736             final Matrix initialMg, final Matrix initialGg, final Matrix accelerometerBias,
1737             final Matrix accelerometerMa, final TurntableGyroscopeCalibratorListener listener) {
1738         this(convertPosition(position), turntableRotationRate, timeInterval, measurements, initialBias, initialMg,
1739                 initialGg, accelerometerBias, accelerometerMa, listener);
1740     }
1741 
1742     /**
1743      * Constructor.
1744      *
1745      * @param position                      position where body kinematics
1746      *                                      measures have been taken.
1747      * @param turntableRotationRate         constant rotation rate at which
1748      *                                      the turntable is spinning. Must
1749      *                                      be expressed in radians per
1750      *                                      second (rad/s).
1751      * @param timeInterval                  time interval between measurements
1752      *                                      being captured expressed in
1753      *                                      seconds (s).
1754      * @param measurements                  collection of body kinematics
1755      *                                      measurements with standard
1756      *                                      deviations taken at the same
1757      *                                      position with zero velocity
1758      *                                      and unknown different
1759      *                                      orientations.
1760      * @param commonAxisUsed                indicates whether z-axis is
1761      *                                      assumed to be common for
1762      *                                      accelerometer and gyroscope.
1763      * @param estimateGDependentCrossBiases true if G-dependent cross biases
1764      *                                      will be estimated, false
1765      *                                      otherwise.
1766      * @param initialBias                   initial gyroscope bias to be
1767      *                                      used to find a solution. This
1768      *                                      must be 3x1 and is expressed in
1769      *                                      radians per second (rad/s).
1770      * @param initialMg                     initial gyroscope scale factors
1771      *                                      and cross coupling errors matrix.
1772      *                                      Must be 3x3.
1773      * @param initialGg                     initial gyroscope G-dependent
1774      *                                      cross biases introduced on the
1775      *                                      gyroscope by the specific
1776      *                                      forces sensed by the
1777      *                                      accelerometer. Must be 3x3.
1778      * @throws IllegalArgumentException if any of the provided values does
1779      *                                  not have proper size or if either
1780      *                                  turntable rotation rate or
1781      *                                  time interval is zero or negative.
1782      */
1783     public TurntableGyroscopeCalibrator(
1784             final NEDPosition position, final double turntableRotationRate, final double timeInterval,
1785             final Collection<StandardDeviationBodyKinematics> measurements, final boolean commonAxisUsed,
1786             final boolean estimateGDependentCrossBiases, final Matrix initialBias, final Matrix initialMg,
1787             final Matrix initialGg) {
1788         this(convertPosition(position), turntableRotationRate, timeInterval, measurements, commonAxisUsed,
1789                 estimateGDependentCrossBiases, initialBias, initialMg, initialGg);
1790     }
1791 
1792     /**
1793      * Constructor.
1794      *
1795      * @param position                      position where body kinematics
1796      *                                      measures have been taken.
1797      * @param turntableRotationRate         constant rotation rate at which
1798      *                                      the turntable is spinning. Must
1799      *                                      be expressed in radians per
1800      *                                      second (rad/s).
1801      * @param timeInterval                  time interval between measurements
1802      *                                      being captured expressed in
1803      *                                      seconds (s).
1804      * @param measurements                  collection of body kinematics
1805      *                                      measurements with standard
1806      *                                      deviations taken at the same
1807      *                                      position with zero velocity and
1808      *                                      unknown different orientations.
1809      * @param commonAxisUsed                indicates whether z-axis is
1810      *                                      assumed to be common for
1811      *                                      accelerometer and gyroscope.
1812      * @param estimateGDependentCrossBiases true if G-dependent cross
1813      *                                      biases will be estimated, false
1814      *                                      otherwise.
1815      * @param initialBias                   initial gyroscope bias to be
1816      *                                      used to find a solution. This
1817      *                                      must be 3x1 and is expressed in
1818      *                                      radians per second (rad/s).
1819      * @param initialMg                     initial gyroscope scale factors
1820      *                                      and cross coupling errors
1821      *                                      matrix. Must be 3x3.
1822      * @param initialGg                     initial gyroscope G-dependent
1823      *                                      cross biases introduced on the
1824      *                                      gyroscope by the specific
1825      *                                      forces sensed by the
1826      *                                      accelerometer. Must be 3x3.
1827      * @param listener                      listener to handle events
1828      *                                      raised by this calibrator.
1829      * @throws IllegalArgumentException if any of the provided values does
1830      *                                  not have proper size or if either
1831      *                                  turntable rotation rate or
1832      *                                  time interval is zero or negative.
1833      */
1834     public TurntableGyroscopeCalibrator(
1835             final NEDPosition position, final double turntableRotationRate, final double timeInterval,
1836             final Collection<StandardDeviationBodyKinematics> measurements, final boolean commonAxisUsed,
1837             final boolean estimateGDependentCrossBiases, final Matrix initialBias, final Matrix initialMg,
1838             final Matrix initialGg, final TurntableGyroscopeCalibratorListener listener) {
1839         this(convertPosition(position), turntableRotationRate, timeInterval, measurements, commonAxisUsed,
1840                 estimateGDependentCrossBiases, initialBias, initialMg, initialGg, listener);
1841     }
1842 
1843     /**
1844      * Constructor.
1845      *
1846      * @param position                      position where body kinematics
1847      *                                      measures have been taken.
1848      * @param turntableRotationRate         constant rotation rate at which
1849      *                                      the turntable is spinning. Must
1850      *                                      be expressed in radians per
1851      *                                      second (rad/s).
1852      * @param timeInterval                  time interval between measurements
1853      *                                      being captured expressed in
1854      *                                      seconds (s).
1855      * @param measurements                  collection of body kinematics
1856      *                                      measurements with standard
1857      *                                      deviations taken at the same
1858      *                                      position with zero velocity
1859      *                                      and unknown different
1860      *                                      orientations.
1861      * @param commonAxisUsed                indicates whether z-axis is
1862      *                                      assumed to be common for
1863      *                                      accelerometer and gyroscope.
1864      * @param estimateGDependentCrossBiases true if G-dependent cross biases
1865      *                                      will be estimated, false
1866      *                                      otherwise.
1867      * @param initialBias                   initial gyroscope bias to be
1868      *                                      used to find a solution. This
1869      *                                      must have length 3 and is
1870      *                                      expressed in radians per second
1871      *                                      (rad/s).
1872      * @param initialMg                     initial gyroscope scale factors
1873      *                                      and cross coupling errors matrix.
1874      *                                      Must be 3x3.
1875      * @param initialGg                     initial gyroscope G-dependent
1876      *                                      cross biases introduced on the
1877      *                                      gyroscope by the specific forces
1878      *                                      sensed by the accelerometer.
1879      *                                      Must be 3x3.
1880      * @throws IllegalArgumentException if any of the provided values does
1881      *                                  not have proper size or if either
1882      *                                  turntable rotation rate or
1883      *                                  time interval is zero or negative.
1884      */
1885     public TurntableGyroscopeCalibrator(
1886             final NEDPosition position, final double turntableRotationRate, final double timeInterval,
1887             final Collection<StandardDeviationBodyKinematics> measurements, final boolean commonAxisUsed,
1888             final boolean estimateGDependentCrossBiases, final double[] initialBias, final Matrix initialMg,
1889             final Matrix initialGg) {
1890         this(convertPosition(position), turntableRotationRate, timeInterval, measurements, commonAxisUsed,
1891                 estimateGDependentCrossBiases, initialBias, initialMg, initialGg);
1892     }
1893 
1894     /**
1895      * Constructor.
1896      *
1897      * @param position                      position where body kinematics
1898      *                                      measures have been taken.
1899      * @param turntableRotationRate         constant rotation rate at which
1900      *                                      the turntable is spinning. Must
1901      *                                      be expressed in radians per
1902      *                                      second (rad/s).
1903      * @param timeInterval                  time interval between measurements
1904      *                                      being captured expressed in
1905      *                                      seconds (s).
1906      * @param measurements                  collection of body kinematics
1907      *                                      measurements with standard
1908      *                                      deviations taken at the same
1909      *                                      position with zero velocity
1910      *                                      and unknown different
1911      *                                      orientations.
1912      * @param commonAxisUsed                indicates whether z-axis is
1913      *                                      assumed to be common for
1914      *                                      accelerometer and gyroscope.
1915      * @param estimateGDependentCrossBiases true if G-dependent cross biases
1916      *                                      will be estimated, false
1917      *                                      otherwise.
1918      * @param initialBias                   initial gyroscope bias to be
1919      *                                      used to find a solution. This
1920      *                                      must have length 3 and is
1921      *                                      expressed in radians per second
1922      *                                      (rad/s).
1923      * @param initialMg                     initial gyroscope scale factors
1924      *                                      and cross coupling errors
1925      *                                      matrix. Must be 3x3.
1926      * @param initialGg                     initial gyroscope G-dependent
1927      *                                      cross biases introduced on the
1928      *                                      gyroscope by the specific forces
1929      *                                      sensed by the accelerometer.
1930      *                                      Must be 3x3.
1931      * @param listener                      listener to handle events raised
1932      *                                      by this calibrator.
1933      * @throws IllegalArgumentException if any of the provided values does
1934      *                                  not have proper size or if either
1935      *                                  turntable rotation rate or
1936      *                                  time interval is zero or negative.
1937      */
1938     public TurntableGyroscopeCalibrator(
1939             final NEDPosition position, final double turntableRotationRate, final double timeInterval,
1940             final Collection<StandardDeviationBodyKinematics> measurements, final boolean commonAxisUsed,
1941             final boolean estimateGDependentCrossBiases, final double[] initialBias, final Matrix initialMg,
1942             final Matrix initialGg, final TurntableGyroscopeCalibratorListener listener) {
1943         this(convertPosition(position), turntableRotationRate, timeInterval, measurements, commonAxisUsed,
1944                 estimateGDependentCrossBiases, initialBias, initialMg, initialGg, listener);
1945     }
1946 
1947     /**
1948      * Constructor.
1949      *
1950      * @param position                      position where body kinematics
1951      *                                      measures have been taken.
1952      * @param turntableRotationRate         constant rotation rate at which
1953      *                                      the turntable is spinning. Must
1954      *                                      be expressed in radians per
1955      *                                      second (rad/s).
1956      * @param timeInterval                  time interval between measurements
1957      *                                      being captured expressed in
1958      *                                      seconds (s).
1959      * @param measurements                  collection of body kinematics
1960      *                                      measurements with standard
1961      *                                      deviations taken at the same
1962      *                                      position with zero velocity
1963      *                                      and unknown different
1964      *                                      orientations.
1965      * @param commonAxisUsed                indicates whether z-axis is
1966      *                                      assumed to be common for
1967      *                                      accelerometer and gyroscope.
1968      * @param estimateGDependentCrossBiases true if G-dependent cross
1969      *                                      biases will be estimated,
1970      *                                      false otherwise.
1971      * @param initialBias                   initial gyroscope bias to be
1972      *                                      used to find a solution. This
1973      *                                      must have length 3 and is
1974      *                                      expressed in radians per second
1975      *                                      (rad/s).
1976      * @param initialMg                     initial gyroscope scale factors
1977      *                                      and cross coupling errors
1978      *                                      matrix. Must be 3x3.
1979      * @param initialGg                     initial gyroscope G-dependent
1980      *                                      cross biases introduced on the
1981      *                                      gyroscope by the specific forces
1982      *                                      sensed by the accelerometer.
1983      *                                      Must be 3x3.
1984      * @param accelerometerBias             known accelerometer bias. This
1985      *                                      must have length 3 and is
1986      *                                      expressed in meters per squared
1987      *                                      second (m/s^2).
1988      * @param accelerometerMa               known accelerometer scale factors
1989      *                                      and cross coupling matrix. Must
1990      *                                      be 3x3.
1991      * @throws IllegalArgumentException if any of the provided values does
1992      *                                  not have proper size or if either
1993      *                                  turntable rotation rate or
1994      *                                  time interval is zero or negative.
1995      */
1996     public TurntableGyroscopeCalibrator(
1997             final NEDPosition position, final double turntableRotationRate, final double timeInterval,
1998             final Collection<StandardDeviationBodyKinematics> measurements, final boolean commonAxisUsed,
1999             final boolean estimateGDependentCrossBiases, final double[] initialBias, final Matrix initialMg,
2000             final Matrix initialGg, final double[] accelerometerBias, final Matrix accelerometerMa) {
2001         this(convertPosition(position), turntableRotationRate, timeInterval, measurements, commonAxisUsed,
2002                 estimateGDependentCrossBiases, initialBias, initialMg, initialGg, accelerometerBias, accelerometerMa);
2003     }
2004 
2005     /**
2006      * Constructor.
2007      *
2008      * @param position                      position where body kinematics
2009      *                                      measures have been taken.
2010      * @param turntableRotationRate         constant rotation rate at which
2011      *                                      the turntable is spinning. Must
2012      *                                      be expressed in radians per
2013      *                                      second (rad/s).
2014      * @param timeInterval                  time interval between measurements
2015      *                                      being captured expressed in
2016      *                                      seconds (s).
2017      * @param measurements                  collection of body kinematics
2018      *                                      measurements with standard
2019      *                                      deviations taken at the same
2020      *                                      position with zero velocity
2021      *                                      and unknown different
2022      *                                      orientations.
2023      * @param commonAxisUsed                indicates whether z-axis is
2024      *                                      assumed to be common for
2025      *                                      accelerometer and gyroscope.
2026      * @param estimateGDependentCrossBiases true if G-dependent cross biases
2027      *                                      will be estimated, false
2028      *                                      otherwise.
2029      * @param initialBias                   initial gyroscope bias to be used
2030      *                                      to find a solution. This must
2031      *                                      have length 3 and is expressed
2032      *                                      in radians per second (rad/s).
2033      * @param initialMg                     initial gyroscope scale factors
2034      *                                      and cross coupling errors matrix.
2035      *                                      Must be 3x3.
2036      * @param initialGg                     initial gyroscope G-dependent
2037      *                                      cross biases introduced on the
2038      *                                      gyroscope by the specific forces
2039      *                                      sensed by the accelerometer. Must
2040      *                                      be 3x3.
2041      * @param accelerometerBias             known accelerometer bias. This
2042      *                                      must have length 3 and is
2043      *                                      expressed in meters per squared
2044      *                                      second (m/s^2).
2045      * @param accelerometerMa               known accelerometer scale factors
2046      *                                      and cross coupling matrix. Must
2047      *                                      be 3x3.
2048      * @param listener                      listener to handle events raised
2049      *                                      by this calibrator.
2050      * @throws IllegalArgumentException if any of the provided values does
2051      *                                  not have proper size or if either
2052      *                                  turntable rotation rate or
2053      *                                  time interval is zero or negative.
2054      */
2055     public TurntableGyroscopeCalibrator(
2056             final NEDPosition position, final double turntableRotationRate, final double timeInterval,
2057             final Collection<StandardDeviationBodyKinematics> measurements, final boolean commonAxisUsed,
2058             final boolean estimateGDependentCrossBiases, final double[] initialBias, final Matrix initialMg,
2059             final Matrix initialGg, final double[] accelerometerBias, final Matrix accelerometerMa,
2060             final TurntableGyroscopeCalibratorListener listener) {
2061         this(convertPosition(position), turntableRotationRate, timeInterval, measurements, commonAxisUsed,
2062                 estimateGDependentCrossBiases, initialBias, initialMg, initialGg, accelerometerBias, accelerometerMa,
2063                 listener);
2064     }
2065 
2066     /**
2067      * Constructor.
2068      *
2069      * @param position                      position where body kinematics
2070      *                                      measures have been taken.
2071      * @param turntableRotationRate         constant rotation rate at which
2072      *                                      the turntable is spinning. Must
2073      *                                      be expressed in radians per
2074      *                                      second (rad/s).
2075      * @param timeInterval                  time interval between measurements
2076      *                                      being captured expressed in
2077      *                                      seconds (s).
2078      * @param measurements                  collection of body kinematics
2079      *                                      measurements with standard
2080      *                                      deviations taken at the same
2081      *                                      position with zero velocity and
2082      *                                      unknown different orientations.
2083      * @param commonAxisUsed                indicates whether z-axis is
2084      *                                      assumed to be common for
2085      *                                      accelerometer and gyroscope.
2086      * @param estimateGDependentCrossBiases true if G-dependent cross biases
2087      *                                      will be estimated, false
2088      *                                      otherwise.
2089      * @param initialBias                   initial gyroscope bias to be
2090      *                                      used to find a solution. This
2091      *                                      must be 3x1 and is expressed in
2092      *                                      radians per second (rad/s).
2093      * @param initialMg                     initial gyroscope scale factors
2094      *                                      and cross coupling errors matrix.
2095      *                                      Must be 3x3.
2096      * @param initialGg                     initial gyroscope G-dependent
2097      *                                      cross biases introduced on the
2098      *                                      gyroscope by the specific forces
2099      *                                      sensed by the accelerometer. Must
2100      *                                      be 3x3.
2101      * @param accelerometerBias             known accelerometer bias. This
2102      *                                      must have length 3 and is
2103      *                                      expressed in meters per squared
2104      *                                      second (m/s^2).
2105      * @param accelerometerMa               known accelerometer scale factors
2106      *                                      and cross coupling matrix. Must
2107      *                                      be 3x3.
2108      * @throws IllegalArgumentException if any of the provided values does
2109      *                                  not have proper size or if either
2110      *                                  turntable rotation rate or
2111      *                                  time interval is zero or negative.
2112      */
2113     public TurntableGyroscopeCalibrator(
2114             final NEDPosition position, final double turntableRotationRate, final double timeInterval,
2115             final Collection<StandardDeviationBodyKinematics> measurements, final boolean commonAxisUsed,
2116             final boolean estimateGDependentCrossBiases, final Matrix initialBias, final Matrix initialMg,
2117             final Matrix initialGg, final Matrix accelerometerBias, final Matrix accelerometerMa) {
2118         this(convertPosition(position), turntableRotationRate, timeInterval, measurements, commonAxisUsed,
2119                 estimateGDependentCrossBiases, initialBias, initialMg, initialGg, accelerometerBias, accelerometerMa);
2120     }
2121 
2122     /**
2123      * Constructor.
2124      *
2125      * @param position                      position where body kinematics
2126      *                                      measures have been taken.
2127      * @param turntableRotationRate         constant rotation rate at which
2128      *                                      the turntable is spinning. Must
2129      *                                      be expressed in radians per
2130      *                                      second (rad/s).
2131      * @param timeInterval                  time interval between measurements
2132      *                                      being captured expressed in
2133      *                                      seconds (s).
2134      * @param measurements                  collection of body kinematics
2135      *                                      measurements with standard
2136      *                                      deviations taken at the same
2137      *                                      position with zero velocity and
2138      *                                      unknown different orientations.
2139      * @param commonAxisUsed                indicates whether z-axis is
2140      *                                      assumed to be common for
2141      *                                      accelerometer and gyroscope.
2142      * @param estimateGDependentCrossBiases true if G-dependent cross biases
2143      *                                      will be estimated, false
2144      *                                      otherwise.
2145      * @param initialBias                   initial gyroscope bias to be used
2146      *                                      to find a solution. This must be
2147      *                                      3x1 and is expressed in radians
2148      *                                      per second (rad/s).
2149      * @param initialMg                     initial gyroscope scale factors
2150      *                                      and cross coupling errors matrix.
2151      *                                      Must be 3x3.
2152      * @param initialGg                     initial gyroscope G-dependent
2153      *                                      cross biases introduced on the
2154      *                                      gyroscope by the specific forces
2155      *                                      sensed by the accelerometer. Must
2156      *                                      be 3x3.
2157      * @param accelerometerBias             known accelerometer bias. This
2158      *                                      must have length 3 and is
2159      *                                      expressed in meters per squared
2160      *                                      second (m/s^2).
2161      * @param accelerometerMa               known accelerometer scale factors
2162      *                                      and cross coupling matrix. Must
2163      *                                      be 3x3.
2164      * @param listener                      listener to handle events raised
2165      *                                      by this calibrator.
2166      * @throws IllegalArgumentException if any of the provided values does
2167      *                                  not have proper size or if either
2168      *                                  turntable rotation rate or
2169      *                                  time interval is zero or negative.
2170      */
2171     public TurntableGyroscopeCalibrator(
2172             final NEDPosition position, final double turntableRotationRate, final double timeInterval,
2173             final Collection<StandardDeviationBodyKinematics> measurements, final boolean commonAxisUsed,
2174             final boolean estimateGDependentCrossBiases, final Matrix initialBias, final Matrix initialMg,
2175             final Matrix initialGg, final Matrix accelerometerBias, final Matrix accelerometerMa,
2176             final TurntableGyroscopeCalibratorListener listener) {
2177         this(convertPosition(position), turntableRotationRate, timeInterval, measurements, commonAxisUsed,
2178                 estimateGDependentCrossBiases, initialBias, initialMg, initialGg, accelerometerBias, accelerometerMa,
2179                 listener);
2180     }
2181 
2182     /**
2183      * Gets known x-coordinate of accelerometer bias to be used to fix
2184      * measured specific force and find cross biases introduced by the
2185      * accelerometer.
2186      * This is expressed in meters per squared second (m/s^2).
2187      *
2188      * @return known x-coordinate of accelerometer bias.
2189      */
2190     @Override
2191     public double getAccelerometerBiasX() {
2192         return accelerometerBiasX;
2193     }
2194 
2195     /**
2196      * Sets known x-coordinate of accelerometer bias to be used to fix
2197      * measured specific force and find cross biases introduced by the
2198      * accelerometer.
2199      * This is expressed in meters per squared second (m/s^2).
2200      *
2201      * @param accelerometerBiasX known x-coordinate of accelerometer bias.
2202      * @throws LockedException if calibrator is currently running.
2203      */
2204     @Override
2205     public void setAccelerometerBiasX(final double accelerometerBiasX) throws LockedException {
2206         if (running) {
2207             throw new LockedException();
2208         }
2209         this.accelerometerBiasX = accelerometerBiasX;
2210     }
2211 
2212     /**
2213      * Gets known y-coordinate of accelerometer bias to be used to fix
2214      * measured specific force and find cross biases introduced by the
2215      * accelerometer.
2216      * This is expressed in meters per squared second (m/s^2).
2217      *
2218      * @return known y-coordinate of accelerometer bias.
2219      */
2220     @Override
2221     public double getAccelerometerBiasY() {
2222         return accelerometerBiasY;
2223     }
2224 
2225     /**
2226      * Sets known y-coordinate of accelerometer bias to be used to fix
2227      * measured specific force and find cross biases introduced by the
2228      * accelerometer.
2229      * This is expressed in meters per squared second (m/s^2).
2230      *
2231      * @param accelerometerBiasY known y-coordinate of accelerometer bias.
2232      * @throws LockedException if calibrator is currently running.
2233      */
2234     @Override
2235     public void setAccelerometerBiasY(final double accelerometerBiasY) throws LockedException {
2236         if (running) {
2237             throw new LockedException();
2238         }
2239         this.accelerometerBiasY = accelerometerBiasY;
2240     }
2241 
2242     /**
2243      * Gets known z-coordinate of accelerometer bias to be used to fix
2244      * measured specific force and find cross biases introduced by the
2245      * accelerometer.
2246      * This is expressed in meters per squared second (m/s^2).
2247      *
2248      * @return known z-coordinate of accelerometer bias.
2249      */
2250     @Override
2251     public double getAccelerometerBiasZ() {
2252         return accelerometerBiasZ;
2253     }
2254 
2255     /**
2256      * Sets known z-coordinate of accelerometer bias to be used to fix
2257      * measured specific force and find cross biases introduced by the
2258      * accelerometer.
2259      * This is expressed in meters per squared second (m/s^2).
2260      *
2261      * @param accelerometerBiasZ known z-coordinate of accelerometer bias.
2262      * @throws LockedException if calibrator is currently running.
2263      */
2264     @Override
2265     public void setAccelerometerBiasZ(final double accelerometerBiasZ) throws LockedException {
2266         if (running) {
2267             throw new LockedException();
2268         }
2269         this.accelerometerBiasZ = accelerometerBiasZ;
2270     }
2271 
2272     /**
2273      * Gets known x-coordinate of accelerometer bias to be used to fix
2274      * measured specific force and find cross biases introduced by the
2275      * accelerometer.
2276      *
2277      * @return known x-coordinate of accelerometer bias.
2278      */
2279     @Override
2280     public Acceleration getAccelerometerBiasXAsAcceleration() {
2281         return new Acceleration(accelerometerBiasX, AccelerationUnit.METERS_PER_SQUARED_SECOND);
2282     }
2283 
2284     /**
2285      * Gets known x-coordinate of accelerometer bias to be used to fix
2286      * measured specific force and find cross biases introduced by the
2287      * accelerometer.
2288      *
2289      * @param result instance where result data will be stored.
2290      */
2291     @Override
2292     public void getAccelerometerBiasXAsAcceleration(final Acceleration result) {
2293         result.setValue(accelerometerBiasX);
2294         result.setUnit(AccelerationUnit.METERS_PER_SQUARED_SECOND);
2295     }
2296 
2297     /**
2298      * Sets known x-coordinate of accelerometer bias to be used to fix
2299      * measured specific force and find cross biases introduced by the
2300      * accelerometer.
2301      *
2302      * @param accelerometerBiasX x-coordinate of accelerometer bias.
2303      * @throws LockedException if calibrator is currently running.
2304      */
2305     @Override
2306     public void setAccelerometerBiasX(final Acceleration accelerometerBiasX) throws LockedException {
2307         if (running) {
2308             throw new LockedException();
2309         }
2310         this.accelerometerBiasX = convertAcceleration(accelerometerBiasX);
2311     }
2312 
2313     /**
2314      * Gets known y-coordinate of accelerometer bias to be used to fix
2315      * measured specific force and find cross biases introduced by the
2316      * accelerometer.
2317      *
2318      * @return known y-coordinate of accelerometer bias.
2319      */
2320     @Override
2321     public Acceleration getAccelerometerBiasYAsAcceleration() {
2322         return new Acceleration(accelerometerBiasY, AccelerationUnit.METERS_PER_SQUARED_SECOND);
2323     }
2324 
2325     /**
2326      * Gets known y-coordinate of accelerometer bias to be used to fix
2327      * measured specific force and find cross biases introduced by the
2328      * accelerometer.
2329      *
2330      * @param result instance where result data will be stored.
2331      */
2332     @Override
2333     public void getAccelerometerBiasYAsAcceleration(final Acceleration result) {
2334         result.setValue(accelerometerBiasY);
2335         result.setUnit(AccelerationUnit.METERS_PER_SQUARED_SECOND);
2336     }
2337 
2338     /**
2339      * Sets known y-coordinate of accelerometer bias to be used to fix
2340      * measured specific force and find cross biases introduced by the
2341      * accelerometer.
2342      *
2343      * @param accelerometerBiasY y-coordinate of accelerometer bias.
2344      * @throws LockedException if calibrator is currently running.
2345      */
2346     @Override
2347     public void setAccelerometerBiasY(final Acceleration accelerometerBiasY) throws LockedException {
2348         if (running) {
2349             throw new LockedException();
2350         }
2351         this.accelerometerBiasY = convertAcceleration(accelerometerBiasY);
2352     }
2353 
2354     /**
2355      * Gets known z-coordinate of accelerometer bias to be used to fix
2356      * measured specific force and find cross biases introduced by the
2357      * accelerometer.
2358      *
2359      * @return known z-coordinate of accelerometer bias.
2360      */
2361     @Override
2362     public Acceleration getAccelerometerBiasZAsAcceleration() {
2363         return new Acceleration(accelerometerBiasZ, AccelerationUnit.METERS_PER_SQUARED_SECOND);
2364     }
2365 
2366     /**
2367      * Gets known z-coordinate of accelerometer bias to be used to fix
2368      * measured specific force and find cross biases introduced by the
2369      * accelerometer.
2370      *
2371      * @param result instance where result data will be stored.
2372      */
2373     @Override
2374     public void getAccelerometerBiasZAsAcceleration(final Acceleration result) {
2375         result.setValue(accelerometerBiasZ);
2376         result.setUnit(AccelerationUnit.METERS_PER_SQUARED_SECOND);
2377     }
2378 
2379     /**
2380      * Sets known z-coordinate of accelerometer bias to be used to fix
2381      * measured specific force and find cross biases introduced by the
2382      * accelerometer.
2383      *
2384      * @param accelerometerBiasZ z-coordinate of accelerometer bias.
2385      * @throws LockedException if calibrator is currently running.
2386      */
2387     @Override
2388     public void setAccelerometerBiasZ(final Acceleration accelerometerBiasZ) throws LockedException {
2389         if (running) {
2390             throw new LockedException();
2391         }
2392         this.accelerometerBiasZ = convertAcceleration(accelerometerBiasZ);
2393     }
2394 
2395     /**
2396      * Sets known accelerometer bias to be used to fix measured specific
2397      * force and find cross biases introduced by the accelerometer.
2398      * This is expressed in meters per squared second (m/s^2).
2399      *
2400      * @param accelerometerBiasX x-coordinate of accelerometer bias.
2401      * @param accelerometerBiasY y-coordinate of accelerometer bias.
2402      * @param accelerometerBiasZ z-coordinate of accelerometer bias.
2403      * @throws LockedException if calibrator is currently running.
2404      */
2405     @Override
2406     public void setAccelerometerBias(
2407             final double accelerometerBiasX, final double accelerometerBiasY, final double accelerometerBiasZ)
2408             throws LockedException {
2409         if (running) {
2410             throw new LockedException();
2411         }
2412 
2413         this.accelerometerBiasX = accelerometerBiasX;
2414         this.accelerometerBiasY = accelerometerBiasY;
2415         this.accelerometerBiasZ = accelerometerBiasZ;
2416     }
2417 
2418     /**
2419      * Sets known accelerometer bias to be used to fix measured specific
2420      * force and find cross biases introduced by the accelerometer.
2421      *
2422      * @param accelerometerBiasX x-coordinate of accelerometer bias.
2423      * @param accelerometerBiasY y-coordinate of accelerometer bias.
2424      * @param accelerometerBiasZ z-coordinate of accelerometer bias.
2425      * @throws LockedException if calibrator is currently running.
2426      */
2427     @Override
2428     public void setAccelerometerBias(
2429             final Acceleration accelerometerBiasX, final Acceleration accelerometerBiasY,
2430             final Acceleration accelerometerBiasZ) throws LockedException {
2431         if (running) {
2432             throw new LockedException();
2433         }
2434 
2435         this.accelerometerBiasX = convertAcceleration(accelerometerBiasX);
2436         this.accelerometerBiasY = convertAcceleration(accelerometerBiasY);
2437         this.accelerometerBiasZ = convertAcceleration(accelerometerBiasZ);
2438     }
2439 
2440     /**
2441      * Gets known accelerometer bias to be used to fix measured specific
2442      * force and find cross biases introduced by the accelerometer.
2443      * This is expressed in meters per squared second (m/s^2).
2444      *
2445      * @return known accelerometer bias.
2446      */
2447     @Override
2448     public double[] getAccelerometerBias() {
2449         final var result = new double[BodyKinematics.COMPONENTS];
2450         getAccelerometerBias(result);
2451         return result;
2452     }
2453 
2454     /**
2455      * Gets known accelerometer bias to be used to fix measured specific
2456      * force and find cross biases introduced by the accelerometer.
2457      * This is expressed in meters per squared second (m/s^2).
2458      *
2459      * @param result instance where result data will be copied to.
2460      * @throws IllegalArgumentException if provided array does not have
2461      *                                  length 3.
2462      */
2463     @Override
2464     public void getAccelerometerBias(final double[] result) {
2465         if (result.length != BodyKinematics.COMPONENTS) {
2466             throw new IllegalArgumentException();
2467         }
2468 
2469         result[0] = accelerometerBiasX;
2470         result[1] = accelerometerBiasY;
2471         result[2] = accelerometerBiasZ;
2472     }
2473 
2474     /**
2475      * Sets known accelerometer bias to be used to fix measured specific
2476      * force and find cross biases introduced by the accelerometer.
2477      * This is expressed in meters per squared second (m/s^2).
2478      *
2479      * @param accelerometerBias known accelerometer bias.
2480      * @throws LockedException          if calibrator is currently running.
2481      * @throws IllegalArgumentException if provided array does not have
2482      *                                  length 3.
2483      */
2484     @Override
2485     public void setAccelerometerBias(final double[] accelerometerBias) throws LockedException {
2486         if (running) {
2487             throw new LockedException();
2488         }
2489 
2490         if (accelerometerBias.length != BodyKinematics.COMPONENTS) {
2491             throw new IllegalArgumentException();
2492         }
2493 
2494         accelerometerBiasX = accelerometerBias[0];
2495         accelerometerBiasY = accelerometerBias[1];
2496         accelerometerBiasZ = accelerometerBias[2];
2497     }
2498 
2499     /**
2500      * Gets known accelerometer bias to be used to fix measured specific
2501      * force and find cross biases introduced by the accelerometer.
2502      * This is expressed in meters per squared second (m/s^2).
2503      *
2504      * @return known accelerometer bias.
2505      */
2506     @Override
2507     public Matrix getAccelerometerBiasAsMatrix() {
2508         Matrix result;
2509         try {
2510             result = new Matrix(BodyKinematics.COMPONENTS, 1);
2511             getAccelerometerBiasAsMatrix(result);
2512         } catch (final WrongSizeException ignore) {
2513             // never happens
2514             result = null;
2515         }
2516         return result;
2517     }
2518 
2519     /**
2520      * Gets known accelerometer bias to be used to fix measured specific
2521      * force and find cross biases introduced by the accelerometer.
2522      * This is expressed in meters per squared second (m/s^2).
2523      *
2524      * @param result instance where result data will be copied to.
2525      * @throws IllegalArgumentException if provided matrix is not 3x1.
2526      */
2527     @Override
2528     public void getAccelerometerBiasAsMatrix(final Matrix result) {
2529         if (result.getRows() != BodyKinematics.COMPONENTS || result.getColumns() != 1) {
2530             throw new IllegalArgumentException();
2531         }
2532         result.setElementAtIndex(0, accelerometerBiasX);
2533         result.setElementAtIndex(1, accelerometerBiasY);
2534         result.setElementAtIndex(2, accelerometerBiasZ);
2535     }
2536 
2537     /**
2538      * Sets known accelerometer bias to be used to fix measured specific
2539      * force and find cross biases introduced by the accelerometer.
2540      * This is expressed in meters per squared second (m/s^2).
2541      *
2542      * @param accelerometerBias known accelerometer bias. Must be 3x1.
2543      * @throws LockedException          if calibrator is currently running.
2544      * @throws IllegalArgumentException if provided matrix is not 3x1.
2545      */
2546     @Override
2547     public void setAccelerometerBias(final Matrix accelerometerBias) throws LockedException {
2548         if (running) {
2549             throw new LockedException();
2550         }
2551         if (accelerometerBias.getRows() != BodyKinematics.COMPONENTS || accelerometerBias.getColumns() != 1) {
2552             throw new IllegalArgumentException();
2553         }
2554 
2555         accelerometerBiasX = accelerometerBias.getElementAtIndex(0);
2556         accelerometerBiasY = accelerometerBias.getElementAtIndex(1);
2557         accelerometerBiasZ = accelerometerBias.getElementAtIndex(2);
2558     }
2559 
2560     /**
2561      * Gets known accelerometer x scaling factor to be used to fix measured
2562      * specific force and find cross biases introduced by the accelerometer.
2563      *
2564      * @return known accelerometer x scaling factor.
2565      */
2566     @Override
2567     public double getAccelerometerSx() {
2568         return accelerometerSx;
2569     }
2570 
2571     /**
2572      * Sets known accelerometer x scaling factor to be used to fix measured
2573      * specific force and find cross biases introduced by the accelerometer.
2574      *
2575      * @param accelerometerSx known accelerometer x scaling factor.
2576      * @throws LockedException if calibrator is currently running.
2577      */
2578     @Override
2579     public void setAccelerometerSx(final double accelerometerSx) throws LockedException {
2580         if (running) {
2581             throw new LockedException();
2582         }
2583         this.accelerometerSx = accelerometerSx;
2584     }
2585 
2586     /**
2587      * Gets known accelerometer y scaling factor to be used to fix measured
2588      * specific force and find cross biases introduced by the accelerometer.
2589      *
2590      * @return known accelerometer y scaling factor.
2591      */
2592     @Override
2593     public double getAccelerometerSy() {
2594         return accelerometerSy;
2595     }
2596 
2597     /**
2598      * Sets known accelerometer y scaling factor to be used to fix measured
2599      * specific force and find cross biases introduced by the accelerometer.
2600      *
2601      * @param accelerometerSy known accelerometer y scaling factor.
2602      * @throws LockedException if calibrator is currently running.
2603      */
2604     @Override
2605     public void setAccelerometerSy(final double accelerometerSy) throws LockedException {
2606         if (running) {
2607             throw new LockedException();
2608         }
2609         this.accelerometerSy = accelerometerSy;
2610     }
2611 
2612     /**
2613      * Gets known accelerometer z scaling factor to be used to fix measured
2614      * specific force and find cross biases introduced by the accelerometer.
2615      *
2616      * @return known accelerometer z scaling factor.
2617      */
2618     @Override
2619     public double getAccelerometerSz() {
2620         return accelerometerSz;
2621     }
2622 
2623     /**
2624      * Sets known accelerometer z scaling factor to be used to fix measured
2625      * specific force and find cross biases introduced by the accelerometer.
2626      *
2627      * @param accelerometerSz known accelerometer z scaling factor.
2628      * @throws LockedException if calibrator is currently running.
2629      */
2630     @Override
2631     public void setAccelerometerSz(final double accelerometerSz) throws LockedException {
2632         if (running) {
2633             throw new LockedException();
2634         }
2635         this.accelerometerSz = accelerometerSz;
2636     }
2637 
2638     /**
2639      * Gets known accelerometer x-y cross coupling error to be used to fix
2640      * measured specific force and find cross biases introduced by the
2641      * accelerometer.
2642      *
2643      * @return known accelerometer x-y cross coupling error.
2644      */
2645     @Override
2646     public double getAccelerometerMxy() {
2647         return accelerometerMxy;
2648     }
2649 
2650     /**
2651      * Sets known accelerometer x-y cross coupling error to be used to fix
2652      * measured specific force and find cross biases introduced by the
2653      * accelerometer.
2654      *
2655      * @param accelerometerMxy known accelerometer x-y cross coupling error.
2656      * @throws LockedException if calibrator is currently running.
2657      */
2658     @Override
2659     public void setAccelerometerMxy(final double accelerometerMxy) throws LockedException {
2660         if (running) {
2661             throw new LockedException();
2662         }
2663         this.accelerometerMxy = accelerometerMxy;
2664     }
2665 
2666     /**
2667      * Gets known accelerometer x-z cross coupling error to be used to fix
2668      * measured specific force and find cross biases introduced by the
2669      * accelerometer.
2670      *
2671      * @return known accelerometer x-z cross coupling error.
2672      */
2673     @Override
2674     public double getAccelerometerMxz() {
2675         return accelerometerMxz;
2676     }
2677 
2678     /**
2679      * Sets known accelerometer x-z cross coupling error to be used to fix
2680      * measured specific force and find cross biases introduced by the
2681      * accelerometer.
2682      *
2683      * @param accelerometerMxz known accelerometer x-z cross coupling error.
2684      * @throws LockedException if calibrator is currently running.
2685      */
2686     @Override
2687     public void setAccelerometerMxz(final double accelerometerMxz) throws LockedException {
2688         if (running) {
2689             throw new LockedException();
2690         }
2691         this.accelerometerMxz = accelerometerMxz;
2692     }
2693 
2694     /**
2695      * Gets known accelerometer y-x cross coupling error to be used to fix
2696      * measured specific force and find cross biases introduced by the
2697      * accelerometer.
2698      *
2699      * @return known accelerometer y-x cross coupling error.
2700      */
2701     @Override
2702     public double getAccelerometerMyx() {
2703         return accelerometerMyx;
2704     }
2705 
2706     /**
2707      * Sets known accelerometer y-x cross coupling error to be used to fix
2708      * measured specific force and find cross biases introduced by the
2709      * accelerometer.
2710      *
2711      * @param accelerometerMyx known accelerometer y-x cross coupling
2712      *                         error.
2713      * @throws LockedException if calibrator is currently running.
2714      */
2715     @Override
2716     public void setAccelerometerMyx(final double accelerometerMyx) throws LockedException {
2717         if (running) {
2718             throw new LockedException();
2719         }
2720         this.accelerometerMyx = accelerometerMyx;
2721     }
2722 
2723     /**
2724      * Gets known accelerometer y-z cross coupling error to be used to fix
2725      * measured specific force and find cross biases introduced by the
2726      * accelerometer.
2727      *
2728      * @return known accelerometer y-z cross coupling error.
2729      */
2730     @Override
2731     public double getAccelerometerMyz() {
2732         return accelerometerMyz;
2733     }
2734 
2735     /**
2736      * Sets known accelerometer y-z cross coupling error to be used to fix
2737      * measured specific force and find cross biases introduced by the
2738      * accelerometer.
2739      *
2740      * @param accelerometerMyz known accelerometer y-z cross coupling
2741      *                         error.
2742      * @throws LockedException if calibrator is currently running.
2743      */
2744     @Override
2745     public void setAccelerometerMyz(final double accelerometerMyz) throws LockedException {
2746         if (running) {
2747             throw new LockedException();
2748         }
2749         this.accelerometerMyz = accelerometerMyz;
2750     }
2751 
2752     /**
2753      * Gets known accelerometer z-x cross coupling error to be used to fix
2754      * measured specific force and find cross biases introduced by the
2755      * accelerometer.
2756      *
2757      * @return known accelerometer z-x cross coupling error.
2758      */
2759     @Override
2760     public double getAccelerometerMzx() {
2761         return accelerometerMzx;
2762     }
2763 
2764     /**
2765      * Sets known accelerometer z-x cross coupling error to be used to fix
2766      * measured specific force and find cross biases introduced by the
2767      * accelerometer.
2768      *
2769      * @param accelerometerMzx known accelerometer z-x cross coupling
2770      *                         error.
2771      * @throws LockedException if calibrator is currently running.
2772      */
2773     @Override
2774     public void setAccelerometerMzx(final double accelerometerMzx) throws LockedException {
2775         if (running) {
2776             throw new LockedException();
2777         }
2778         this.accelerometerMzx = accelerometerMzx;
2779     }
2780 
2781     /**
2782      * Gets known accelerometer z-y cross coupling error to be used to fix
2783      * measured specific force and find cross biases introduced by the
2784      * accelerometer.
2785      *
2786      * @return known accelerometer z-y cross coupling error.
2787      */
2788     @Override
2789     public double getAccelerometerMzy() {
2790         return accelerometerMzy;
2791     }
2792 
2793     /**
2794      * Sets known accelerometer z-y cross coupling error to be used to fix
2795      * measured specific force and find cross biases introduced by the
2796      * accelerometer.
2797      *
2798      * @param accelerometerMzy known accelerometer z-y cross coupling
2799      *                         error.
2800      * @throws LockedException if calibrator is currently running.
2801      */
2802     @Override
2803     public void setAccelerometerMzy(final double accelerometerMzy) throws LockedException {
2804         if (running) {
2805             throw new LockedException();
2806         }
2807         this.accelerometerMzy = accelerometerMzy;
2808     }
2809 
2810     /**
2811      * Sets known accelerometer scaling factors to be used to fix measured
2812      * specific force and find cross biases introduced by the
2813      * accelerometer.
2814      *
2815      * @param accelerometerSx known accelerometer x scaling factor.
2816      * @param accelerometerSy known accelerometer y scaling factor.
2817      * @param accelerometerSz known accelerometer z scaling factor.
2818      * @throws LockedException if calibrator is currently running.
2819      */
2820     @Override
2821     public void setAccelerometerScalingFactors(
2822             final double accelerometerSx, final double accelerometerSy, final double accelerometerSz)
2823             throws LockedException {
2824         if (running) {
2825             throw new LockedException();
2826         }
2827         this.accelerometerSx = accelerometerSx;
2828         this.accelerometerSy = accelerometerSy;
2829         this.accelerometerSz = accelerometerSz;
2830     }
2831 
2832     /**
2833      * Sets known accelerometer cross coupling errors to be used to fix
2834      * measured specific force and find cross biases introduced by the
2835      * accelerometer.
2836      *
2837      * @param accelerometerMxy known accelerometer x-y cross coupling
2838      *                         error.
2839      * @param accelerometerMxz known accelerometer x-z cross coupling
2840      *                         error.
2841      * @param accelerometerMyx known accelerometer y-x cross coupling
2842      *                         error.
2843      * @param accelerometerMyz known accelerometer y-z cross coupling
2844      *                         error.
2845      * @param accelerometerMzx known accelerometer z-x cross coupling
2846      *                         error.
2847      * @param accelerometerMzy known accelerometer z-y cross coupling
2848      *                         error.
2849      * @throws LockedException if calibrator is currently running.
2850      */
2851     @Override
2852     public void setAccelerometerCrossCouplingErrors(
2853             final double accelerometerMxy, final double accelerometerMxz, final double accelerometerMyx,
2854             final double accelerometerMyz, final double accelerometerMzx, final double accelerometerMzy)
2855             throws LockedException {
2856         if (running) {
2857             throw new LockedException();
2858         }
2859         this.accelerometerMxy = accelerometerMxy;
2860         this.accelerometerMxz = accelerometerMxz;
2861         this.accelerometerMyx = accelerometerMyx;
2862         this.accelerometerMyz = accelerometerMyz;
2863         this.accelerometerMzx = accelerometerMzx;
2864         this.accelerometerMzy = accelerometerMzy;
2865     }
2866 
2867     /**
2868      * Sets known accelerometer scaling factors and cross coupling errors
2869      * to be used to fix measured specific force and find cross biases
2870      * introduced by the accelerometer.
2871      *
2872      * @param accelerometerSx  known accelerometer x scaling factor.
2873      * @param accelerometerSy  known accelerometer y scaling factor.
2874      * @param accelerometerSz  known accelerometer z scaling factor.
2875      * @param accelerometerMxy known accelerometer x-y cross coupling
2876      *                         error.
2877      * @param accelerometerMxz known accelerometer x-z cross coupling
2878      *                         error.
2879      * @param accelerometerMyx known accelerometer y-x cross coupling
2880      *                         error.
2881      * @param accelerometerMyz known accelerometer y-z cross coupling
2882      *                         error.
2883      * @param accelerometerMzx known accelerometer z-x cross coupling
2884      *                         error.
2885      * @param accelerometerMzy known accelerometer z-y cross coupling
2886      *                         error.
2887      * @throws LockedException if calibrator is currently running.
2888      */
2889     @Override
2890     public void setAccelerometerScalingFactorsAndCrossCouplingErrors(
2891             final double accelerometerSx, final double accelerometerSy, final double accelerometerSz,
2892             final double accelerometerMxy, final double accelerometerMxz, final double accelerometerMyx,
2893             final double accelerometerMyz, final double accelerometerMzx, final double accelerometerMzy)
2894             throws LockedException {
2895         if (running) {
2896             throw new LockedException();
2897         }
2898         setAccelerometerScalingFactors(accelerometerSx, accelerometerSy, accelerometerSz);
2899         setAccelerometerCrossCouplingErrors(accelerometerMxy, accelerometerMxz, accelerometerMyx,
2900                 accelerometerMyz, accelerometerMzx, accelerometerMzy);
2901     }
2902 
2903     /**
2904      * Gets known accelerometer scale factors and cross coupling
2905      * errors matrix.
2906      *
2907      * @return known accelerometer scale factors and cross coupling
2908      * errors matrix.
2909      */
2910     @Override
2911     public Matrix getAccelerometerMa() {
2912         Matrix result;
2913         try {
2914             result = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
2915             getAccelerometerMa(result);
2916         } catch (final WrongSizeException ignore) {
2917             // never happens
2918             result = null;
2919         }
2920         return result;
2921     }
2922 
2923     /**
2924      * Gets known accelerometer scale factors and cross coupling
2925      * errors matrix.
2926      *
2927      * @param result instance where data will be stored.
2928      * @throws IllegalArgumentException if provided matrix is not 3x3.
2929      */
2930     @Override
2931     public void getAccelerometerMa(final Matrix result) {
2932         if (result.getRows() != BodyKinematics.COMPONENTS || result.getColumns() != BodyKinematics.COMPONENTS) {
2933             throw new IllegalArgumentException();
2934         }
2935         result.setElementAtIndex(0, accelerometerSx);
2936         result.setElementAtIndex(1, accelerometerMyx);
2937         result.setElementAtIndex(2, accelerometerMzx);
2938 
2939         result.setElementAtIndex(3, accelerometerMxy);
2940         result.setElementAtIndex(4, accelerometerSy);
2941         result.setElementAtIndex(5, accelerometerMzy);
2942 
2943         result.setElementAtIndex(6, accelerometerMxz);
2944         result.setElementAtIndex(7, accelerometerMyz);
2945         result.setElementAtIndex(8, accelerometerSz);
2946     }
2947 
2948     /**
2949      * Sets known accelerometer scale factors and cross coupling
2950      * errors matrix.
2951      *
2952      * @param accelerometerMa known accelerometer scale factors and
2953      *                        cross coupling errors matrix. Must be 3x3.
2954      * @throws LockedException          if calibrator is currently running.
2955      * @throws IllegalArgumentException if provided matrix is not 3x3.
2956      */
2957     @Override
2958     public void setAccelerometerMa(final Matrix accelerometerMa) throws LockedException {
2959         if (running) {
2960             throw new LockedException();
2961         }
2962         if (accelerometerMa.getRows() != BodyKinematics.COMPONENTS
2963                 || accelerometerMa.getColumns() != BodyKinematics.COMPONENTS) {
2964             throw new IllegalArgumentException();
2965         }
2966 
2967         accelerometerSx = accelerometerMa.getElementAtIndex(0);
2968         accelerometerMyx = accelerometerMa.getElementAtIndex(1);
2969         accelerometerMzx = accelerometerMa.getElementAtIndex(2);
2970 
2971         accelerometerMxy = accelerometerMa.getElementAtIndex(3);
2972         accelerometerSy = accelerometerMa.getElementAtIndex(4);
2973         accelerometerMzy = accelerometerMa.getElementAtIndex(5);
2974 
2975         accelerometerMxz = accelerometerMa.getElementAtIndex(6);
2976         accelerometerMyz = accelerometerMa.getElementAtIndex(7);
2977         accelerometerSz = accelerometerMa.getElementAtIndex(8);
2978     }
2979 
2980     /**
2981      * Gets initial x-coordinate of gyroscope bias to be used to find
2982      * a solution.
2983      * This is expressed in radians per second (rad/s).
2984      *
2985      * @return initial x-coordinate of gyroscope bias.
2986      */
2987     public double getInitialBiasX() {
2988         return initialBiasX;
2989     }
2990 
2991     /**
2992      * Sets initial x-coordinate of gyroscope bias to be used to find
2993      * a solution.
2994      * This is expressed in radians per second (rad/s).
2995      *
2996      * @param initialBiasX initial x-coordinate of gyroscope bias.
2997      * @throws LockedException if calibrator is currently running.
2998      */
2999     public void setInitialBiasX(final double initialBiasX) throws LockedException {
3000         if (running) {
3001             throw new LockedException();
3002         }
3003         this.initialBiasX = initialBiasX;
3004     }
3005 
3006     /**
3007      * Gets initial y-coordinate of gyroscope bias to be used to find
3008      * a solution.
3009      * This is expressed in radians per second (rad/s).
3010      *
3011      * @return initial y-coordinate of gyroscope bias.
3012      */
3013     public double getInitialBiasY() {
3014         return initialBiasY;
3015     }
3016 
3017     /**
3018      * Sets initial y-coordinate of gyroscope bias to be used to find
3019      * a solution.
3020      * This is expressed in radians per second (rad/s).
3021      *
3022      * @param initialBiasY initial y-coordinate of gyroscope bias.
3023      * @throws LockedException if calibrator is currently running.
3024      */
3025     public void setInitialBiasY(final double initialBiasY) throws LockedException {
3026         if (running) {
3027             throw new LockedException();
3028         }
3029         this.initialBiasY = initialBiasY;
3030     }
3031 
3032     /**
3033      * Gets initial z-coordinate of gyroscope bias ot be used to find
3034      * a solution.
3035      * This is expressed in radians per second (rad/s).
3036      *
3037      * @return initial z-coordinate of gyroscope bias.
3038      */
3039     public double getInitialBiasZ() {
3040         return initialBiasZ;
3041     }
3042 
3043     /**
3044      * Sets initial z-coordinate of gyroscope bias to be used to find
3045      * a solution.
3046      * This is expressed in radians per second (rad/s).
3047      *
3048      * @param initialBiasZ initial z-coordinate of gyroscope bias.
3049      * @throws LockedException if calibrator is currently running.
3050      */
3051     public void setInitialBiasZ(final double initialBiasZ) throws LockedException {
3052         if (running) {
3053             throw new LockedException();
3054         }
3055         this.initialBiasZ = initialBiasZ;
3056     }
3057 
3058     /**
3059      * Gets initial x-coordinate of gyroscope bias to be used to find a
3060      * solution.
3061      *
3062      * @return initial x-coordinate of gyroscope bias.
3063      */
3064     public AngularSpeed getInitialBiasAngularSpeedX() {
3065         return new AngularSpeed(initialBiasX, AngularSpeedUnit.RADIANS_PER_SECOND);
3066     }
3067 
3068     /**
3069      * Gets initial x-coordinate of gyroscope bias to be used to find a
3070      * solution.
3071      *
3072      * @param result instance where result data will be stored.
3073      */
3074     public void getInitialBiasAngularSpeedX(final AngularSpeed result) {
3075         result.setValue(initialBiasX);
3076         result.setUnit(AngularSpeedUnit.RADIANS_PER_SECOND);
3077     }
3078 
3079     /**
3080      * Sets initial x-coordinate of gyroscope bias to be used to find a
3081      * solution.
3082      *
3083      * @param initialBiasX initial x-coordinate of gyroscope bias.
3084      * @throws LockedException if calibrator is currently running.
3085      */
3086     public void setInitialBiasX(final AngularSpeed initialBiasX) throws LockedException {
3087         if (running) {
3088             throw new LockedException();
3089         }
3090         this.initialBiasX = convertAngularSpeed(initialBiasX);
3091     }
3092 
3093     /**
3094      * Gets initial y-coordinate of gyroscope bias to be used to find a
3095      * solution.
3096      *
3097      * @return initial y-coordinate of gyroscope bias.
3098      */
3099     public AngularSpeed getInitialBiasAngularSpeedY() {
3100         return new AngularSpeed(initialBiasY, AngularSpeedUnit.RADIANS_PER_SECOND);
3101     }
3102 
3103     /**
3104      * Gets initial y-coordinate of gyroscope bias to be used to find a
3105      * solution.
3106      *
3107      * @param result instance where result data will be stored.
3108      */
3109     public void getInitialBiasAngularSpeedY(final AngularSpeed result) {
3110         result.setValue(initialBiasY);
3111         result.setUnit(AngularSpeedUnit.RADIANS_PER_SECOND);
3112     }
3113 
3114     /**
3115      * Sets initial y-coordinate of gyroscope bias to be used to find a
3116      * solution.
3117      *
3118      * @param initialBiasY initial y-coordinate of gyroscope bias.
3119      * @throws LockedException if calibrator is currently running.
3120      */
3121     public void setInitialBiasY(final AngularSpeed initialBiasY) throws LockedException {
3122         if (running) {
3123             throw new LockedException();
3124         }
3125         this.initialBiasY = convertAngularSpeed(initialBiasY);
3126     }
3127 
3128     /**
3129      * Gets initial z-coordinate of gyroscope bias to be used to find a
3130      * solution.
3131      *
3132      * @return initial z-coordinate of gyroscope bias.
3133      */
3134     public AngularSpeed getInitialBiasAngularSpeedZ() {
3135         return new AngularSpeed(initialBiasZ, AngularSpeedUnit.RADIANS_PER_SECOND);
3136     }
3137 
3138     /**
3139      * Gets initial z-coordinate of gyroscope bias to be used to find a
3140      * solution.
3141      *
3142      * @param result instance where result data will be stored.
3143      */
3144     public void getInitialBiasAngularSpeedZ(final AngularSpeed result) {
3145         result.setValue(initialBiasZ);
3146         result.setUnit(AngularSpeedUnit.RADIANS_PER_SECOND);
3147     }
3148 
3149     /**
3150      * Sets initial z-coordinate of gyroscope bias to be used to find a
3151      * solution.
3152      *
3153      * @param initialBiasZ initial z-coordinate of gyroscope bias.
3154      * @throws LockedException if calibrator is currently running.
3155      */
3156     public void setInitialBiasZ(final AngularSpeed initialBiasZ) throws LockedException {
3157         if (running) {
3158             throw new LockedException();
3159         }
3160         this.initialBiasZ = convertAngularSpeed(initialBiasZ);
3161     }
3162 
3163     /**
3164      * Sets initial bias coordinates of gyroscope used to find a solution
3165      * expressed in radians per second (rad/s).
3166      *
3167      * @param initialBiasX initial x-coordinate of gyroscope bias.
3168      * @param initialBiasY initial y-coordinate of gyroscope bias.
3169      * @param initialBiasZ initial z-coordinate of gyroscope bias.
3170      * @throws LockedException if calibrator is currently running.
3171      */
3172     public void setInitialBias(
3173             final double initialBiasX, final double initialBiasY, final double initialBiasZ) throws LockedException {
3174         if (running) {
3175             throw new LockedException();
3176         }
3177         this.initialBiasX = initialBiasX;
3178         this.initialBiasY = initialBiasY;
3179         this.initialBiasZ = initialBiasZ;
3180     }
3181 
3182     /**
3183      * Sets initial bias coordinates of gyroscope used to find a solution.
3184      *
3185      * @param initialBiasX initial x-coordinate of gyroscope bias.
3186      * @param initialBiasY initial y-coordinate of gyroscope bias.
3187      * @param initialBiasZ initial z-coordinate of gyroscope bias.
3188      * @throws LockedException if calibrator is currently running.
3189      */
3190     public void setInitialBias(
3191             final AngularSpeed initialBiasX, final AngularSpeed initialBiasY, final AngularSpeed initialBiasZ)
3192             throws LockedException {
3193         if (running) {
3194             throw new LockedException();
3195         }
3196         this.initialBiasX = convertAngularSpeed(initialBiasX);
3197         this.initialBiasY = convertAngularSpeed(initialBiasY);
3198         this.initialBiasZ = convertAngularSpeed(initialBiasZ);
3199     }
3200 
3201     /**
3202      * Gets initial x scaling factor of gyroscope.
3203      *
3204      * @return initial x scaling factor of gyroscope.
3205      */
3206     @Override
3207     public double getInitialSx() {
3208         return initialSx;
3209     }
3210 
3211     /**
3212      * Sets initial x scaling factor of gyroscope.
3213      *
3214      * @param initialSx initial x scaling factor of gyroscope.
3215      * @throws LockedException if calibrator is currently running.
3216      */
3217     @Override
3218     public void setInitialSx(final double initialSx) throws LockedException {
3219         if (running) {
3220             throw new LockedException();
3221         }
3222         this.initialSx = initialSx;
3223     }
3224 
3225     /**
3226      * Gets initial y scaling factor of gyroscope.
3227      *
3228      * @return initial y scaling factor of gyroscope.
3229      */
3230     @Override
3231     public double getInitialSy() {
3232         return initialSy;
3233     }
3234 
3235     /**
3236      * Sets initial y scaling factor of gyroscope.
3237      *
3238      * @param initialSy initial y scaling factor of gyroscope.
3239      * @throws LockedException if calibrator is currently running.
3240      */
3241     @Override
3242     public void setInitialSy(final double initialSy) throws LockedException {
3243         if (running) {
3244             throw new LockedException();
3245         }
3246         this.initialSy = initialSy;
3247     }
3248 
3249     /**
3250      * Gets initial z scaling factor of gyroscope.
3251      *
3252      * @return initial z scaling factor of gyroscope.
3253      */
3254     @Override
3255     public double getInitialSz() {
3256         return initialSz;
3257     }
3258 
3259     /**
3260      * Sets initial z scaling factor of gyroscope.
3261      *
3262      * @param initialSz initial z scaling factor of gyroscope.
3263      * @throws LockedException if calibrator is currently running.
3264      */
3265     @Override
3266     public void setInitialSz(final double initialSz) throws LockedException {
3267         if (running) {
3268             throw new LockedException();
3269         }
3270         this.initialSz = initialSz;
3271     }
3272 
3273     /**
3274      * Gets initial x-y cross coupling error of gyroscope.
3275      *
3276      * @return initial x-y cross coupling error of gyroscope.
3277      */
3278     @Override
3279     public double getInitialMxy() {
3280         return initialMxy;
3281     }
3282 
3283     /**
3284      * Sets initial x-y cross coupling error of gyroscope.
3285      *
3286      * @param initialMxy initial x-y cross coupling error of gyroscope.
3287      * @throws LockedException if calibrator is currently running.
3288      */
3289     @Override
3290     public void setInitialMxy(final double initialMxy) throws LockedException {
3291         if (running) {
3292             throw new LockedException();
3293         }
3294         this.initialMxy = initialMxy;
3295     }
3296 
3297     /**
3298      * Gets initial x-z cross coupling error of gyroscope.
3299      *
3300      * @return initial x-z cross coupling error of gyroscope.
3301      */
3302     @Override
3303     public double getInitialMxz() {
3304         return initialMxz;
3305     }
3306 
3307     /**
3308      * Sets initial x-z cross coupling error of gyroscope.
3309      *
3310      * @param initialMxz initial x-z cross coupling error of gyroscope.
3311      * @throws LockedException if calibrator is currently running.
3312      */
3313     @Override
3314     public void setInitialMxz(final double initialMxz) throws LockedException {
3315         if (running) {
3316             throw new LockedException();
3317         }
3318         this.initialMxz = initialMxz;
3319     }
3320 
3321     /**
3322      * Gets initial y-x cross coupling error of gyroscope.
3323      *
3324      * @return initial y-x cross coupling error of gyroscope.
3325      */
3326     @Override
3327     public double getInitialMyx() {
3328         return initialMyx;
3329     }
3330 
3331     /**
3332      * Sets initial y-x cross coupling error of gyroscope.
3333      *
3334      * @param initialMyx initial y-x cross coupling error of gyroscope.
3335      * @throws LockedException if calibrator is currently running.
3336      */
3337     @Override
3338     public void setInitialMyx(final double initialMyx) throws LockedException {
3339         if (running) {
3340             throw new LockedException();
3341         }
3342         this.initialMyx = initialMyx;
3343     }
3344 
3345     /**
3346      * Gets initial y-z cross coupling error of gyroscope.
3347      *
3348      * @return initial y-z cross coupling error of gyroscope.
3349      */
3350     @Override
3351     public double getInitialMyz() {
3352         return initialMyz;
3353     }
3354 
3355     /**
3356      * Sets initial y-z cross coupling error of gyroscope.
3357      *
3358      * @param initialMyz initial y-z cross coupling error of gyroscope.
3359      * @throws LockedException if calibrator is currently running.
3360      */
3361     @Override
3362     public void setInitialMyz(final double initialMyz) throws LockedException {
3363         if (running) {
3364             throw new LockedException();
3365         }
3366         this.initialMyz = initialMyz;
3367     }
3368 
3369     /**
3370      * Gets initial z-x cross coupling error of gyroscope.
3371      *
3372      * @return initial z-x cross coupling error of gyroscope.
3373      */
3374     @Override
3375     public double getInitialMzx() {
3376         return initialMzx;
3377     }
3378 
3379     /**
3380      * Sets initial z-x cross coupling error of gyroscope.
3381      *
3382      * @param initialMzx initial z-x cross coupling error of gyroscope.
3383      * @throws LockedException if calibrator is currently running.
3384      */
3385     @Override
3386     public void setInitialMzx(final double initialMzx) throws LockedException {
3387         if (running) {
3388             throw new LockedException();
3389         }
3390         this.initialMzx = initialMzx;
3391     }
3392 
3393     /**
3394      * Gets initial z-y cross coupling error of gyroscope.
3395      *
3396      * @return initial z-y cross coupling error of gyroscope.
3397      */
3398     @Override
3399     public double getInitialMzy() {
3400         return initialMzy;
3401     }
3402 
3403     /**
3404      * Sets initial z-y cross coupling error of gyroscope.
3405      *
3406      * @param initialMzy initial z-y cross coupling error of gyroscope.
3407      * @throws LockedException if calibrator is currently running.
3408      */
3409     @Override
3410     public void setInitialMzy(final double initialMzy) throws LockedException {
3411         if (running) {
3412             throw new LockedException();
3413         }
3414         this.initialMzy = initialMzy;
3415     }
3416 
3417     /**
3418      * Sets initial scaling factors of gyroscope.
3419      *
3420      * @param initialSx initial x scaling factor of gyroscope.
3421      * @param initialSy initial y scaling factor of gyroscope.
3422      * @param initialSz initial z scaling factor of gyroscope.
3423      * @throws LockedException if calibrator is currently running.
3424      */
3425     @Override
3426     public void setInitialScalingFactors(
3427             final double initialSx, final double initialSy, final double initialSz) throws LockedException {
3428         if (running) {
3429             throw new LockedException();
3430         }
3431         this.initialSx = initialSx;
3432         this.initialSy = initialSy;
3433         this.initialSz = initialSz;
3434     }
3435 
3436     /**
3437      * Sets initial cross coupling errors of gyroscope.
3438      *
3439      * @param initialMxy initial x-y cross coupling error of gyroscope.
3440      * @param initialMxz initial x-z cross coupling error of gyroscope.
3441      * @param initialMyx initial y-x cross coupling error of gyroscope.
3442      * @param initialMyz initial y-z cross coupling error of gyroscope.
3443      * @param initialMzx initial z-x cross coupling error of gyroscope.
3444      * @param initialMzy initial z-y cross coupling error of gyroscope.
3445      * @throws LockedException if calibrator is currently running.
3446      */
3447     @Override
3448     public void setInitialCrossCouplingErrors(
3449             final double initialMxy, final double initialMxz, final double initialMyx,
3450             final double initialMyz, final double initialMzx, final double initialMzy) throws LockedException {
3451         if (running) {
3452             throw new LockedException();
3453         }
3454         this.initialMxy = initialMxy;
3455         this.initialMxz = initialMxz;
3456         this.initialMyx = initialMyx;
3457         this.initialMyz = initialMyz;
3458         this.initialMzx = initialMzx;
3459         this.initialMzy = initialMzy;
3460     }
3461 
3462     /**
3463      * Sets initial scaling factors and cross coupling errors of
3464      * gyroscope.
3465      *
3466      * @param initialSx  initial x scaling factor of gyroscope.
3467      * @param initialSy  initial y scaling factor of gyroscope.
3468      * @param initialSz  initial z scaling factor of gyroscope.
3469      * @param initialMxy initial x-y cross coupling error of gyroscope.
3470      * @param initialMxz initial x-z cross coupling error of gyroscope.
3471      * @param initialMyx initial y-x cross coupling error of gyroscope.
3472      * @param initialMyz initial y-z cross coupling error of gyroscope.
3473      * @param initialMzx initial z-x cross coupling error of gyroscope.
3474      * @param initialMzy initial z-y cross coupling error of gyroscope.
3475      * @throws LockedException if calibrator is currently running.
3476      */
3477     @Override
3478     public void setInitialScalingFactorsAndCrossCouplingErrors(
3479             final double initialSx, final double initialSy, final double initialSz,
3480             final double initialMxy, final double initialMxz, final double initialMyx,
3481             final double initialMyz, final double initialMzx, final double initialMzy) throws LockedException {
3482         if (running) {
3483             throw new LockedException();
3484         }
3485         setInitialScalingFactors(initialSx, initialSy, initialSz);
3486         setInitialCrossCouplingErrors(initialMxy, initialMxz, initialMyx, initialMyz, initialMzx, initialMzy);
3487     }
3488 
3489     /**
3490      * Gets initial gyroscope bias to be used to find a solution as
3491      * an array.
3492      * Array values are expressed in radians per second (rad/s).
3493      *
3494      * @return array containing coordinates of initial gyroscope bias.
3495      */
3496     public double[] getInitialBias() {
3497         final var result = new double[BodyKinematics.COMPONENTS];
3498         getInitialBias(result);
3499         return result;
3500     }
3501 
3502     /**
3503      * Gets initial gyroscope bias to be used to find a solution as
3504      * an array.
3505      * Array values are expressed in radians per second (rad/s).
3506      *
3507      * @param result instance where result data will be copied to.
3508      * @throws IllegalArgumentException if provided array does not have length 3.
3509      */
3510     public void getInitialBias(final double[] result) {
3511         if (result.length != BodyKinematics.COMPONENTS) {
3512             throw new IllegalArgumentException();
3513         }
3514         result[0] = initialBiasX;
3515         result[1] = initialBiasY;
3516         result[2] = initialBiasZ;
3517     }
3518 
3519     /**
3520      * Sets initial gyroscope bias to be used to find a solution as
3521      * an array.
3522      * Array values are expressed in radians per second (rad/s).
3523      *
3524      * @param initialBias initial bias to find a solution.
3525      * @throws LockedException          if calibrator is currently running.
3526      * @throws IllegalArgumentException if provided array does not have length 3.
3527      */
3528     public void setInitialBias(final double[] initialBias) throws LockedException {
3529         if (running) {
3530             throw new LockedException();
3531         }
3532 
3533         if (initialBias.length != BodyKinematics.COMPONENTS) {
3534             throw new IllegalArgumentException();
3535         }
3536         initialBiasX = initialBias[0];
3537         initialBiasY = initialBias[1];
3538         initialBiasZ = initialBias[2];
3539     }
3540 
3541     /**
3542      * Gets initial gyroscope bias to be used to find a solution as a
3543      * column matrix.
3544      *
3545      * @return initial gyroscope bias to be used to find a solution as a
3546      * column matrix.
3547      */
3548     public Matrix getInitialBiasAsMatrix() {
3549         Matrix result;
3550         try {
3551             result = new Matrix(BodyKinematics.COMPONENTS, 1);
3552             getInitialBiasAsMatrix(result);
3553         } catch (final WrongSizeException ignore) {
3554             // never happens
3555             result = null;
3556         }
3557         return result;
3558     }
3559 
3560     /**
3561      * Gets initial gyroscope bias to be used to find a solution as a
3562      * column matrix.
3563      *
3564      * @param result instance where result data will be copied to.
3565      * @throws IllegalArgumentException if provided matrix is not 3x1.
3566      */
3567     public void getInitialBiasAsMatrix(final Matrix result) {
3568         if (result.getRows() != BodyKinematics.COMPONENTS || result.getColumns() != 1) {
3569             throw new IllegalArgumentException();
3570         }
3571         result.setElementAtIndex(0, initialBiasX);
3572         result.setElementAtIndex(1, initialBiasY);
3573         result.setElementAtIndex(2, initialBiasZ);
3574     }
3575 
3576     /**
3577      * Sets initial gyroscope bias to be used to find a solution as
3578      * an array.
3579      *
3580      * @param initialBias initial gyroscope bias to find a solution.
3581      * @throws LockedException          if calibrator is currently running.
3582      * @throws IllegalArgumentException if provided matrix is not 3x1.
3583      */
3584     public void setInitialBias(final Matrix initialBias) throws LockedException {
3585         if (running) {
3586             throw new LockedException();
3587         }
3588         if (initialBias.getRows() != BodyKinematics.COMPONENTS || initialBias.getColumns() != 1) {
3589             throw new IllegalArgumentException();
3590         }
3591 
3592         initialBiasX = initialBias.getElementAtIndex(0);
3593         initialBiasY = initialBias.getElementAtIndex(1);
3594         initialBiasZ = initialBias.getElementAtIndex(2);
3595     }
3596 
3597     /**
3598      * Gets initial bias coordinates of gyroscope used to find a solution.
3599      *
3600      * @return initial bias coordinates.
3601      */
3602     public AngularSpeedTriad getInitialBiasAsTriad() {
3603         return new AngularSpeedTriad(AngularSpeedUnit.RADIANS_PER_SECOND, initialBiasX, initialBiasY, initialBiasZ);
3604     }
3605 
3606     /**
3607      * Gets initial bias coordinates of gyroscope used to find a solution.
3608      *
3609      * @param result instance where result will be stored.
3610      */
3611     public void getInitialBiasAsTriad(final AngularSpeedTriad result) {
3612         result.setValueCoordinatesAndUnit(initialBiasX, initialBiasY, initialBiasZ,
3613                 AngularSpeedUnit.RADIANS_PER_SECOND);
3614     }
3615 
3616     /**
3617      * Sets initial bias coordinates of gyroscope used to find a solution.
3618      *
3619      * @param initialBias initial bias coordinates to be set.
3620      * @throws LockedException if calibrator is currently running.
3621      */
3622     public void setInitialBias(final AngularSpeedTriad initialBias) throws LockedException {
3623         if (running) {
3624             throw new LockedException();
3625         }
3626 
3627         initialBiasX = convertAngularSpeed(initialBias.getValueX(), initialBias.getUnit());
3628         initialBiasY = convertAngularSpeed(initialBias.getValueY(), initialBias.getUnit());
3629         initialBiasZ = convertAngularSpeed(initialBias.getValueZ(), initialBias.getUnit());
3630     }
3631 
3632     /**
3633      * Gets initial gyroscope scale factors and cross coupling errors
3634      * matrix.
3635      *
3636      * @return initial gyroscope scale factors and cross coupling errors
3637      * matrix.
3638      */
3639     @Override
3640     public Matrix getInitialMg() {
3641         Matrix result;
3642         try {
3643             result = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
3644             getInitialMg(result);
3645         } catch (final WrongSizeException ignore) {
3646             // never happens
3647             result = null;
3648         }
3649         return result;
3650     }
3651 
3652     /**
3653      * Gets initial gyroscope scale factors and cross coupling errors
3654      * matrix.
3655      *
3656      * @param result instance where data will be stored.
3657      * @throws IllegalArgumentException if provided matrix is not 3x3.
3658      */
3659     @Override
3660     public void getInitialMg(final Matrix result) {
3661         if (result.getRows() != BodyKinematics.COMPONENTS || result.getColumns() != BodyKinematics.COMPONENTS) {
3662             throw new IllegalArgumentException();
3663         }
3664         result.setElementAtIndex(0, initialSx);
3665         result.setElementAtIndex(1, initialMyx);
3666         result.setElementAtIndex(2, initialMzx);
3667 
3668         result.setElementAtIndex(3, initialMxy);
3669         result.setElementAtIndex(4, initialSy);
3670         result.setElementAtIndex(5, initialMzy);
3671 
3672         result.setElementAtIndex(6, initialMxz);
3673         result.setElementAtIndex(7, initialMyz);
3674         result.setElementAtIndex(8, initialSz);
3675     }
3676 
3677     /**
3678      * Sets initial gyroscope scale factors and cross coupling errors matrix.
3679      *
3680      * @param initialMg initial scale factors and cross coupling errors matrix.
3681      * @throws IllegalArgumentException if provided matrix is not 3x3.
3682      * @throws LockedException          if calibrator is currently running.
3683      */
3684     @Override
3685     public void setInitialMg(final Matrix initialMg) throws LockedException {
3686         if (running) {
3687             throw new LockedException();
3688         }
3689         if (initialMg.getRows() != BodyKinematics.COMPONENTS || initialMg.getColumns() != BodyKinematics.COMPONENTS) {
3690             throw new IllegalArgumentException();
3691         }
3692 
3693         initialSx = initialMg.getElementAtIndex(0);
3694         initialMyx = initialMg.getElementAtIndex(1);
3695         initialMzx = initialMg.getElementAtIndex(2);
3696 
3697         initialMxy = initialMg.getElementAtIndex(3);
3698         initialSy = initialMg.getElementAtIndex(4);
3699         initialMzy = initialMg.getElementAtIndex(5);
3700 
3701         initialMxz = initialMg.getElementAtIndex(6);
3702         initialMyz = initialMg.getElementAtIndex(7);
3703         initialSz = initialMg.getElementAtIndex(8);
3704     }
3705 
3706     /**
3707      * Gets initial G-dependent cross biases introduced on the gyroscope by the
3708      * specific forces sensed by the accelerometer.
3709      *
3710      * @return a 3x3 matrix containing initial g-dependent cross biases.
3711      */
3712     @Override
3713     public Matrix getInitialGg() {
3714         return new Matrix(initialGg);
3715     }
3716 
3717     /**
3718      * Gets initial G-dependent cross biases introduced on the gyroscope by the
3719      * specific forces sensed by the accelerometer.
3720      *
3721      * @param result instance where data will be stored.
3722      * @throws IllegalArgumentException if provided matrix is not 3x3.
3723      */
3724     @Override
3725     public void getInitialGg(final Matrix result) {
3726         if (result.getRows() != BodyKinematics.COMPONENTS || result.getColumns() != BodyKinematics.COMPONENTS) {
3727             throw new IllegalArgumentException();
3728         }
3729 
3730         result.copyFrom(initialGg);
3731     }
3732 
3733     /**
3734      * Sets initial G-dependent cross biases introduced on the gyroscope by the
3735      * specific forces sensed by the accelerometer.
3736      *
3737      * @param initialGg g-dependent cross biases.
3738      * @throws LockedException          if calibrator is currently running.
3739      * @throws IllegalArgumentException if provided matrix is not 3x3.
3740      */
3741     @Override
3742     public void setInitialGg(final Matrix initialGg) throws LockedException {
3743         if (running) {
3744             throw new LockedException();
3745         }
3746 
3747         if (initialGg.getRows() != BodyKinematics.COMPONENTS || initialGg.getColumns() != BodyKinematics.COMPONENTS) {
3748             throw new IllegalArgumentException();
3749         }
3750 
3751         initialGg.copyTo(this.initialGg);
3752     }
3753 
3754     /**
3755      * Gets constant rotation rate at which the turntable is spinning.
3756      * This is expressed in radians per second (rad/s).
3757      *
3758      * @return constant rotation rate of turntable.
3759      */
3760     public double getTurntableRotationRate() {
3761         return turntableRotationRate;
3762     }
3763 
3764     /**
3765      * Sets constant rotation rate at which the turntable is spinning.
3766      * This is expressed in radians per second (rad/s).
3767      *
3768      * @param turntableRotationRate constant rotation rate of turntable.
3769      * @throws LockedException          if calibrator is currently running
3770      * @throws IllegalArgumentException if provided value is zero or
3771      *                                  negative.
3772      */
3773     public void setTurntableRotationRate(final double turntableRotationRate) throws LockedException {
3774         if (running) {
3775             throw new LockedException();
3776         }
3777         if (turntableRotationRate <= 0.0) {
3778             throw new IllegalArgumentException();
3779         }
3780 
3781         this.turntableRotationRate = turntableRotationRate;
3782     }
3783 
3784     /**
3785      * Gets constant rotation rate at which the turntable is spinning.
3786      *
3787      * @return constant rotation rate of turntable.
3788      */
3789     public AngularSpeed getTurntableRotationRateAsAngularSpeed() {
3790         return new AngularSpeed(turntableRotationRate, AngularSpeedUnit.RADIANS_PER_SECOND);
3791     }
3792 
3793     /**
3794      * Gets constant rotation rate at which the turntable is spinning.
3795      *
3796      * @param result instance where result will be stored.
3797      */
3798     public void getTurntableRotationRateAsAngularSpeed(final AngularSpeed result) {
3799         result.setValue(turntableRotationRate);
3800         result.setUnit(AngularSpeedUnit.RADIANS_PER_SECOND);
3801     }
3802 
3803     /**
3804      * Sets constant rotation rate at which the turntable is spinning.
3805      *
3806      * @param turntableRotationRate constant rotation rate of turntable.
3807      * @throws LockedException          if calibrator is currently running.
3808      * @throws IllegalArgumentException if provided value is zero or
3809      *                                  negative.
3810      */
3811     public void setTurntableRotationRate(final AngularSpeed turntableRotationRate) throws LockedException {
3812         if (running) {
3813             throw new LockedException();
3814         }
3815         setTurntableRotationRate(convertAngularSpeed(turntableRotationRate));
3816     }
3817 
3818     /**
3819      * Gets time interval between measurements being captured expressed in
3820      * seconds (s).
3821      *
3822      * @return time interval between measurements.
3823      */
3824     public double getTimeInterval() {
3825         return timeInterval;
3826     }
3827 
3828     /**
3829      * Sets time interval between measurements being captured expressed in
3830      * seconds (s).
3831      *
3832      * @param timeInterval time interval between measurements.
3833      * @throws LockedException          if calibrator is currently running.
3834      * @throws IllegalArgumentException if provided value is zero or
3835      *                                  negative.
3836      */
3837     public void setTimeInterval(final double timeInterval) throws LockedException {
3838         if (running) {
3839             throw new LockedException();
3840         }
3841 
3842         if (timeInterval <= 0.0) {
3843             throw new IllegalArgumentException();
3844         }
3845         this.timeInterval = timeInterval;
3846     }
3847 
3848     /**
3849      * Gets time interval between measurements being captured.
3850      *
3851      * @return time interval between measurements.
3852      */
3853     public Time getTimeIntervalAsTime() {
3854         return new Time(timeInterval, TimeUnit.SECOND);
3855     }
3856 
3857     /**
3858      * Gets time interval between measurements being captured.
3859      *
3860      * @param result instance where result will be stored.
3861      */
3862     public void getTimeIntervalAsTime(final Time result) {
3863         result.setValue(timeInterval);
3864         result.setUnit(TimeUnit.SECOND);
3865     }
3866 
3867     /**
3868      * Sets time interval between measurements being captured.
3869      *
3870      * @param timeInterval time interval between measurements.
3871      * @throws LockedException if calibrator is currently running.
3872      */
3873     public void setTimeInterval(final Time timeInterval) throws LockedException {
3874         if (running) {
3875             throw new LockedException();
3876         }
3877         setTimeInterval(convertTime(timeInterval));
3878     }
3879 
3880     /**
3881      * Gets a collection of body kinematics measurements taken at
3882      * a given position with different unknown orientations and containing
3883      * the standard deviations of accelerometer and gyroscope measurements.
3884      *
3885      * @return collection of body kinematics measurements at a known position
3886      * with unknown orientations.
3887      */
3888     @Override
3889     public Collection<StandardDeviationBodyKinematics> getMeasurements() {
3890         return measurements;
3891     }
3892 
3893     /**
3894      * Sets a collection of body kinematics measurements taken at
3895      * a given position with different unknown orientations and containing
3896      * the standard deviations of accelerometer and gyroscope measurements.
3897      *
3898      * @param measurements collection of body kinematics measurements at a
3899      *                     known position with unknown orientations.
3900      * @throws LockedException if calibrator is currently running.
3901      */
3902     @Override
3903     public void setMeasurements(final Collection<StandardDeviationBodyKinematics> measurements) throws LockedException {
3904         if (running) {
3905             throw new LockedException();
3906         }
3907         this.measurements = measurements;
3908     }
3909 
3910     /**
3911      * Gets position where body kinematics measures have been taken expressed in
3912      * ECEF coordinates.
3913      *
3914      * @return position where body kinematics measures have been taken.
3915      */
3916     public ECEFPosition getEcefPosition() {
3917         return position;
3918     }
3919 
3920     /**
3921      * Gets position where body kinematics measures have been taken expressed in
3922      * ECEF coordinates.
3923      *
3924      * @param position position where body kinematics measures have been taken.
3925      * @throws LockedException if calibrator is currently running.
3926      */
3927     public void setPosition(final ECEFPosition position) throws LockedException {
3928         if (running) {
3929             throw new LockedException();
3930         }
3931 
3932         this.position = position;
3933     }
3934 
3935     /**
3936      * Gets position where body kinematics measures have been taken expressed in
3937      * NED coordinates.
3938      *
3939      * @return position where body kinematics measures have been taken or null if
3940      * not available.
3941      */
3942     public NEDPosition getNedPosition() {
3943         final var result = new NEDPosition();
3944         return getNedPosition(result) ? result : null;
3945     }
3946 
3947     /**
3948      * Gets position where body kinematics measures have been taken expressed in
3949      * NED coordinates.
3950      *
3951      * @param result instance where result will be stored.
3952      * @return true if NED position could be computed, false otherwise.
3953      */
3954     public boolean getNedPosition(final NEDPosition result) {
3955         if (position != null) {
3956             final var velocity = new NEDVelocity();
3957             ECEFtoNEDPositionVelocityConverter.convertECEFtoNED(
3958                     position.getX(), position.getY(), position.getZ(),
3959                     0.0, 0.0, 0.0, result, velocity);
3960             return true;
3961         } else {
3962             return false;
3963         }
3964     }
3965 
3966     /**
3967      * Sets position where body kinematics measures have been taken expressed in
3968      * NED coordinates.
3969      *
3970      * @param position position where body kinematics measures have been taken.
3971      * @throws LockedException if calibrator is currently running.
3972      */
3973     public void setPosition(final NEDPosition position) throws LockedException {
3974         if (running) {
3975             throw new LockedException();
3976         }
3977 
3978         this.position = convertPosition(position);
3979     }
3980 
3981     /**
3982      * Indicates the type of measurement or sequence used by this calibrator.
3983      *
3984      * @return type of measurement or sequence used by this calibrator.
3985      */
3986     @Override
3987     public GyroscopeCalibratorMeasurementOrSequenceType getMeasurementOrSequenceType() {
3988         return GyroscopeCalibratorMeasurementOrSequenceType.STANDARD_DEVIATION_BODY_KINEMATICS_MEASUREMENT;
3989     }
3990 
3991     /**
3992      * Indicates whether this calibrator requires ordered measurements or sequences
3993      * in a list or not.
3994      *
3995      * @return true if measurements or sequences must be ordered, false otherwise.
3996      */
3997     @Override
3998     public boolean isOrderedMeasurementsOrSequencesRequired() {
3999         return false;
4000     }
4001 
4002     /**
4003      * Indicates whether this calibrator requires quality scores for each
4004      * measurement/sequence or not.
4005      *
4006      * @return true if quality scores are required, false otherwise.
4007      */
4008     @Override
4009     public boolean isQualityScoresRequired() {
4010         return false;
4011     }
4012 
4013     /**
4014      * Indicates whether z-axis is assumed to be common for accelerometer and
4015      * gyroscope.
4016      * When enabled, this eliminates 3 variables from Ma matrix.
4017      *
4018      * @return true if z-axis is assumed to be common for accelerometer and gyroscope,
4019      * false otherwise.
4020      */
4021     @Override
4022     public boolean isCommonAxisUsed() {
4023         return commonAxisUsed;
4024     }
4025 
4026     /**
4027      * Specifies whether z-axis is assumed to be common for accelerometer and
4028      * gyroscope.
4029      * When enabled, this eliminates 3 variables from Ma matrix.
4030      *
4031      * @param commonAxisUsed true if z-axis is assumed to be common for accelerometer
4032      *                       and gyroscope, false otherwise.
4033      * @throws LockedException if calibrator is currently running.
4034      */
4035     @Override
4036     public void setCommonAxisUsed(final boolean commonAxisUsed) throws LockedException {
4037         if (running) {
4038             throw new LockedException();
4039         }
4040 
4041         this.commonAxisUsed = commonAxisUsed;
4042     }
4043 
4044     /**
4045      * Indicates whether G-dependent cross biases are being estimated
4046      * or not.
4047      * When enabled, this adds 9 variables from Gg matrix.
4048      *
4049      * @return true if G-dependent cross biases will be estimated,
4050      * false otherwise.
4051      */
4052     public boolean isGDependentCrossBiasesEstimated() {
4053         return estimateGDependentCrossBiases;
4054     }
4055 
4056     /**
4057      * Specifies whether G-dependent cross biases are being estimated
4058      * or not.
4059      * When enabled, this adds 9 variables from Gg matrix.
4060      *
4061      * @param estimateGDependentCrossBiases true if G-dependent cross
4062      *                                      biases will be estimated,
4063      *                                      false otherwise.
4064      * @throws LockedException if calibrator is currently running.
4065      */
4066     public void setGDependentCrossBiasesEstimated(final boolean estimateGDependentCrossBiases) throws LockedException {
4067         if (running) {
4068             throw new LockedException();
4069         }
4070 
4071         this.estimateGDependentCrossBiases = estimateGDependentCrossBiases;
4072     }
4073 
4074     /**
4075      * Gets listener to handle events raised by this estimator.
4076      *
4077      * @return listener to handle events raised by this estimator.
4078      */
4079     public TurntableGyroscopeCalibratorListener getListener() {
4080         return listener;
4081     }
4082 
4083     /**
4084      * Sets listener to handle events raised by this estimator.
4085      *
4086      * @param listener listener to handle events raised by this estimator.
4087      * @throws LockedException if calibrator is currently running.
4088      */
4089     public void setListener(final TurntableGyroscopeCalibratorListener listener) throws LockedException {
4090         if (running) {
4091             throw new LockedException();
4092         }
4093 
4094         this.listener = listener;
4095     }
4096 
4097     /**
4098      * Gets minimum number of required measurements.
4099      *
4100      * @return minimum number of required measurements.
4101      */
4102     @Override
4103     public int getMinimumRequiredMeasurementsOrSequences() {
4104         if (commonAxisUsed) {
4105             if (estimateGDependentCrossBiases) {
4106                 return MINIMUM_MEASUREMENTS_COMMON_Z_AXIS_AND_CROSS_BIASES;
4107             } else {
4108                 return MINIMUM_MEASUREMENTS_COMMON_Z_AXIS;
4109             }
4110         } else {
4111             if (estimateGDependentCrossBiases) {
4112                 return MINIMUM_MEASUREMENTS_GENERAL_AND_CROSS_BIASES;
4113             } else {
4114                 return MINIMUM_MEASUREMENTS_GENERAL;
4115             }
4116         }
4117     }
4118 
4119     /**
4120      * Indicates whether calibrator is ready to start.
4121      *
4122      * @return true if calibrator is ready, false otherwise.
4123      */
4124     @Override
4125     public boolean isReady() {
4126         return measurements != null && measurements.size() >= getMinimumRequiredMeasurementsOrSequences();
4127     }
4128 
4129     /**
4130      * Indicates whether calibrator is currently running or not.
4131      *
4132      * @return true if calibrator is running, false otherwise.
4133      */
4134     @Override
4135     public boolean isRunning() {
4136         return running;
4137     }
4138 
4139     /**
4140      * Estimates gyroscope calibration parameters containing bias, scale factors,
4141      * cross-coupling errors and G-dependent coupling.
4142      *
4143      * @throws LockedException      if calibrator is currently running.
4144      * @throws NotReadyException    if calibrator is not ready.
4145      * @throws CalibrationException if estimation fails for numerical reasons.
4146      */
4147     @Override
4148     public void calibrate() throws LockedException, NotReadyException, CalibrationException {
4149         if (running) {
4150             throw new LockedException();
4151         }
4152 
4153         if (!isReady()) {
4154             throw new NotReadyException();
4155         }
4156 
4157         try {
4158             running = true;
4159 
4160             if (listener != null) {
4161                 listener.onCalibrateStart(this);
4162             }
4163 
4164             if (commonAxisUsed) {
4165                 if (estimateGDependentCrossBiases) {
4166                     calibrateCommonAxisAndGDependentCrossBiases();
4167                 } else {
4168                     calibrateCommonAxis();
4169                 }
4170             } else {
4171                 if (estimateGDependentCrossBiases) {
4172                     calibrateGeneralAndGDependentCrossBiases();
4173                 } else {
4174                     calibrateGeneral();
4175                 }
4176             }
4177 
4178             if (listener != null) {
4179                 listener.onCalibrateEnd(this);
4180             }
4181 
4182         } catch (final AlgebraException | FittingException | com.irurueta.numerical.NotReadyException |
4183                        InvalidSourceAndDestinationFrameTypeException e) {
4184             throw new CalibrationException(e);
4185         } finally {
4186             running = false;
4187         }
4188     }
4189 
4190     /**
4191      * Gets array containing x,y,z components of estimated gyroscope biases
4192      * expressed in radians per second (rad/s).
4193      *
4194      * @return array containing x,y,z components of estimated gyroscope biases.
4195      */
4196     @Override
4197     public double[] getEstimatedBiases() {
4198         return estimatedBiases;
4199     }
4200 
4201     /**
4202      * Gets array containing x,y,z components of estimated gyroscope biases
4203      * expressed in radians per second (rad/s).
4204      *
4205      * @param result instance where estimated gyroscope biases will be stored.
4206      * @return true if result instance was updated, false otherwise (when estimation
4207      * is not yet available).
4208      */
4209     @Override
4210     public boolean getEstimatedBiases(final double[] result) {
4211         if (estimatedBiases != null) {
4212             System.arraycopy(estimatedBiases, 0, result, 0, estimatedBiases.length);
4213             return true;
4214         } else {
4215             return false;
4216         }
4217     }
4218 
4219     /**
4220      * Gets column matrix containing x,y,z components of estimated gyroscope biases
4221      * expressed in radians per second (rad/s).
4222      *
4223      * @return column matrix containing x,y,z components of estimated gyroscope
4224      * biases.
4225      */
4226     @Override
4227     public Matrix getEstimatedBiasesAsMatrix() {
4228         return estimatedBiases != null ? Matrix.newFromArray(estimatedBiases) : null;
4229     }
4230 
4231     /**
4232      * Gets column matrix containing x,y,z components of estimated gyroscope biases
4233      * expressed in radians per second (rad/s).
4234      *
4235      * @param result instance where result data will be stored.
4236      * @return true if result was updated, false otherwise.
4237      * @throws WrongSizeException if provided result instance has invalid size.
4238      */
4239     @Override
4240     public boolean getEstimatedBiasesAsMatrix(final Matrix result) throws WrongSizeException {
4241         if (estimatedBiases != null) {
4242             result.fromArray(estimatedBiases);
4243             return true;
4244         } else {
4245             return false;
4246         }
4247     }
4248 
4249     /**
4250      * Gets x coordinate of estimated gyroscope bias expressed in radians per
4251      * second (rad/s).
4252      *
4253      * @return x coordinate of estimated gyroscope bias or null if not available.
4254      */
4255     @Override
4256     public Double getEstimatedBiasX() {
4257         return estimatedBiases != null ? estimatedBiases[0] : null;
4258     }
4259 
4260     /**
4261      * Gets y coordinate of estimated gyroscope bias expressed in radians per
4262      * second (rad/s).
4263      *
4264      * @return y coordinate of estimated gyroscope bias or null if not available.
4265      */
4266     @Override
4267     public Double getEstimatedBiasY() {
4268         return estimatedBiases != null ? estimatedBiases[1] : null;
4269     }
4270 
4271     /**
4272      * Gets z coordinate of estimated gyroscope bias expressed in radians per
4273      * second (rad/s).
4274      *
4275      * @return z coordinate of estimated gyroscope bias or null if not available.
4276      */
4277     @Override
4278     public Double getEstimatedBiasZ() {
4279         return estimatedBiases != null ? estimatedBiases[2] : null;
4280     }
4281 
4282     /**
4283      * Gets x coordinate of estimated gyroscope bias.
4284      *
4285      * @return x coordinate of estimated gyroscope bias or null if not available.
4286      */
4287     @Override
4288     public AngularSpeed getEstimatedBiasAngularSpeedX() {
4289         return estimatedBiases != null
4290                 ? new AngularSpeed(estimatedBiases[0], AngularSpeedUnit.RADIANS_PER_SECOND) : null;
4291     }
4292 
4293     /**
4294      * Gets x coordinate of estimated gyroscope bias.
4295      *
4296      * @param result instance where result will be stored.
4297      * @return true if result was updated, false if estimation is not available.
4298      */
4299     @Override
4300     public boolean getEstimatedBiasAngularSpeedX(final AngularSpeed result) {
4301         if (estimatedBiases != null) {
4302             result.setValue(estimatedBiases[0]);
4303             result.setUnit(AngularSpeedUnit.RADIANS_PER_SECOND);
4304             return true;
4305         } else {
4306             return false;
4307         }
4308     }
4309 
4310     /**
4311      * Gets y coordinate of estimated gyroscope bias.
4312      *
4313      * @return y coordinate of estimated gyroscope bias or null if not available.
4314      */
4315     @Override
4316     public AngularSpeed getEstimatedBiasAngularSpeedY() {
4317         return estimatedBiases != null
4318                 ? new AngularSpeed(estimatedBiases[1], AngularSpeedUnit.RADIANS_PER_SECOND) : null;
4319     }
4320 
4321     /**
4322      * Gets y coordinate of estimated gyroscope bias.
4323      *
4324      * @param result instance where result will be stored.
4325      * @return true if result was updated, false if estimation is not available.
4326      */
4327     @Override
4328     public boolean getEstimatedBiasAngularSpeedY(final AngularSpeed result) {
4329         if (estimatedBiases != null) {
4330             result.setValue(estimatedBiases[1]);
4331             result.setUnit(AngularSpeedUnit.RADIANS_PER_SECOND);
4332             return true;
4333         } else {
4334             return false;
4335         }
4336     }
4337 
4338     /**
4339      * Gets z coordinate of estimated gyroscope bias.
4340      *
4341      * @return z coordinate of estimated gyroscope bias or null if not available.
4342      */
4343     @Override
4344     public AngularSpeed getEstimatedBiasAngularSpeedZ() {
4345         return estimatedBiases != null
4346                 ? new AngularSpeed(estimatedBiases[2], AngularSpeedUnit.RADIANS_PER_SECOND) : null;
4347     }
4348 
4349     /**
4350      * Gets z coordinate of estimated gyroscope bias.
4351      *
4352      * @param result instance where result will be stored.
4353      * @return true if result was updated, false if estimation is not available.
4354      */
4355     @Override
4356     public boolean getEstimatedBiasAngularSpeedZ(final AngularSpeed result) {
4357         if (estimatedBiases != null) {
4358             result.setValue(estimatedBiases[2]);
4359             result.setUnit(AngularSpeedUnit.RADIANS_PER_SECOND);
4360             return true;
4361         } else {
4362             return false;
4363         }
4364     }
4365 
4366     /**
4367      * Gets estimated gyroscope bias.
4368      *
4369      * @return estimated gyroscope bias or null if not available.
4370      */
4371     @Override
4372     public AngularSpeedTriad getEstimatedBiasAsTriad() {
4373         return estimatedBiases != null
4374                 ? new AngularSpeedTriad(AngularSpeedUnit.RADIANS_PER_SECOND,
4375                 estimatedBiases[0], estimatedBiases[1], estimatedBiases[2])
4376                 : null;
4377     }
4378 
4379     /**
4380      * Gets estimated gyroscope bias.
4381      *
4382      * @param result instance where result will be stored.
4383      * @return true if estimated gyroscope bias is available and result was
4384      * modified, false otherwise.
4385      */
4386     @Override
4387     public boolean getEstimatedBiasAsTriad(final AngularSpeedTriad result) {
4388         if (estimatedBiases != null) {
4389             result.setValueCoordinatesAndUnit(estimatedBiases[0], estimatedBiases[1], estimatedBiases[2],
4390                     AngularSpeedUnit.RADIANS_PER_SECOND);
4391             return true;
4392         } else {
4393             return false;
4394         }
4395     }
4396 
4397     /**
4398      * Gets estimated gyroscope scale factors and cross coupling errors.
4399      * This is the product of matrix Tg containing cross coupling errors and Kg
4400      * containing scaling factors.
4401      * So that:
4402      * <pre>
4403      *     Mg = [sx    mxy  mxz] = Tg*Kg
4404      *          [myx   sy   myz]
4405      *          [mzx   mzy  sz ]
4406      * </pre>
4407      * Where:
4408      * <pre>
4409      *     Kg = [sx 0   0 ]
4410      *          [0  sy  0 ]
4411      *          [0  0   sz]
4412      * </pre>
4413      * and
4414      * <pre>
4415      *     Tg = [1          -alphaXy    alphaXz ]
4416      *          [alphaYx    1           -alphaYz]
4417      *          [-alphaZx   alphaZy     1       ]
4418      * </pre>
4419      * Hence:
4420      * <pre>
4421      *     Mg = [sx    mxy  mxz] = Tg*Kg =  [sx             -sy * alphaXy   sz * alphaXz ]
4422      *          [myx   sy   myz]            [sx * alphaYx   sy              -sz * alphaYz]
4423      *          [mzx   mzy  sz ]            [-sx * alphaZx  sy * alphaZy    sz           ]
4424      * </pre>
4425      * This instance allows any 3x3 matrix however, typically alphaYx, alphaZx and alphaZy
4426      * are considered to be zero if the gyroscope z-axis is assumed to be the same
4427      * as the body z-axis. When this is assumed, myx = mzx = mzy = 0 and the Mg matrix
4428      * becomes upper diagonal:
4429      * <pre>
4430      *     Mg = [sx    mxy  mxz]
4431      *          [0     sy   myz]
4432      *          [0     0    sz ]
4433      * </pre>
4434      * Values of this matrix are unit-less.
4435      *
4436      * @return estimated gyroscope scale factors and cross coupling errors, or null
4437      * if not available.
4438      */
4439     @Override
4440     public Matrix getEstimatedMg() {
4441         return estimatedMg;
4442     }
4443 
4444     /**
4445      * Gets estimated gyroscope x-axis scale factor.
4446      *
4447      * @return estimated gyroscope x-axis scale factor or null
4448      * if not available.
4449      */
4450     @Override
4451     public Double getEstimatedSx() {
4452         return estimatedMg != null ? estimatedMg.getElementAt(0, 0) : null;
4453     }
4454 
4455     /**
4456      * Gets estimated gyroscope y-axis scale factor.
4457      *
4458      * @return estimated gyroscope y-axis scale factor or null
4459      * if not available.
4460      */
4461     @Override
4462     public Double getEstimatedSy() {
4463         return estimatedMg != null ? estimatedMg.getElementAt(1, 1) : null;
4464     }
4465 
4466     /**
4467      * Gets estimated gyroscope z-axis scale factor.
4468      *
4469      * @return estimated gyroscope z-axis scale factor or null
4470      * if not available.
4471      */
4472     @Override
4473     public Double getEstimatedSz() {
4474         return estimatedMg != null ? estimatedMg.getElementAt(2, 2) : null;
4475     }
4476 
4477     /**
4478      * Gets estimated gyroscope x-y cross-coupling error.
4479      *
4480      * @return estimated gyroscope x-y cross-coupling error or null
4481      * if not available.
4482      */
4483     @Override
4484     public Double getEstimatedMxy() {
4485         return estimatedMg != null ? estimatedMg.getElementAt(0, 1) : null;
4486     }
4487 
4488     /**
4489      * Gets estimated gyroscope x-z cross-coupling error.
4490      *
4491      * @return estimated gyroscope x-z cross-coupling error or null
4492      * if not available.
4493      */
4494     @Override
4495     public Double getEstimatedMxz() {
4496         return estimatedMg != null ? estimatedMg.getElementAt(0, 2) : null;
4497     }
4498 
4499     /**
4500      * Gets estimated gyroscope y-x cross-coupling error.
4501      *
4502      * @return estimated gyroscope y-x cross-coupling error or null
4503      * if not available.
4504      */
4505     @Override
4506     public Double getEstimatedMyx() {
4507         return estimatedMg != null ? estimatedMg.getElementAt(1, 0) : null;
4508     }
4509 
4510     /**
4511      * Gets estimated gyroscope y-z cross-coupling error.
4512      *
4513      * @return estimated gyroscope y-z cross-coupling error or null
4514      * if not available.
4515      */
4516     @Override
4517     public Double getEstimatedMyz() {
4518         return estimatedMg != null ? estimatedMg.getElementAt(1, 2) : null;
4519     }
4520 
4521     /**
4522      * Gets estimated gyroscope z-x cross-coupling error.
4523      *
4524      * @return estimated gyroscope z-x cross-coupling error or null
4525      * if not available.
4526      */
4527     @Override
4528     public Double getEstimatedMzx() {
4529         return estimatedMg != null ? estimatedMg.getElementAt(2, 0) : null;
4530     }
4531 
4532     /**
4533      * Gets estimated gyroscope z-y cross-coupling error.
4534      *
4535      * @return estimated gyroscope z-y cross-coupling error or null
4536      * if not available.
4537      */
4538     @Override
4539     public Double getEstimatedMzy() {
4540         return estimatedMg != null ? estimatedMg.getElementAt(2, 1) : null;
4541     }
4542 
4543     /**
4544      * Gets estimated G-dependent cross biases introduced on the gyroscope by the
4545      * specific forces sensed by the accelerometer.
4546      * This instance allows any 3x3 matrix.
4547      *
4548      * @return estimated G-dependent cross biases.
4549      */
4550     @Override
4551     public Matrix getEstimatedGg() {
4552         return estimatedGg;
4553     }
4554 
4555     /**
4556      * Gets estimated covariance matrix for estimated parameters.
4557      * Diagonal elements of the matrix contains variance for the following
4558      * parameters (following indicated order): bgx, bgy, bgz, sx, sy, sz,
4559      * mxy, mxz, myx, myz, mzx, mzy, gg11, gg21, gg31, gg12, gg22, gg32,
4560      * gg13, gg23, gg33.
4561      *
4562      * @return estimated covariance matrix for estimated parameters.
4563      */
4564     @Override
4565     public Matrix getEstimatedCovariance() {
4566         return estimatedCovariance;
4567     }
4568 
4569     /**
4570      * Gets estimated chi square value.
4571      *
4572      * @return estimated chi square value.
4573      */
4574     @Override
4575     public double getEstimatedChiSq() {
4576         return estimatedChiSq;
4577     }
4578 
4579     /**
4580      * Gets estimated chi square degrees of freedom. Degrees of freedom is equal to the number of sampled data minus the
4581      * number of estimated parameters.
4582      *
4583      * @return estimated degrees of freedom of chi square value
4584      */
4585     @Override
4586     public int getEstimatedChiSqDegreesOfFreedom() {
4587         return estimatedChiSqDegreesOfFreedom;
4588     }
4589 
4590     /**
4591      * Gets estimated reduced chi square value. This is equal to estimated chi square value divided by its degrees of
4592      * freedom. Ideally this value should be close to 1.0, indicating that fit is optimal.
4593      * A value larger than 1.0 indicates that fit is not good or noise has been underestimated, and a value smaller than
4594      * 1.0 indicates that there is overfitting or noise has been overestimated.
4595      *
4596      * @return estimated reduced chi square value
4597      */
4598     @Override
4599     public double getEstimatedReducedChiSq() {
4600         return estimatedReducedChiSq;
4601     }
4602 
4603     /**
4604      * Gets estimated mean square error respect to provided measurements.
4605      *
4606      * @return estimated mean square error respect to provided measurements.
4607      */
4608     @Override
4609     public double getEstimatedMse() {
4610         return estimatedMse;
4611     }
4612 
4613     /**
4614      * Gets estimated probability of finding a smaller chi square value expressed as a value between 0.0 and 1.0. The
4615      * smaller the found chi square value is, the better the fit of the estimated parameters to the actual parameter.
4616      * Thus, the smaller the chance of finding a smaller chi square value, then the better the estimated fit is.
4617      *
4618      * @return estimated probability of finding a smaller chi square value.
4619      */
4620     @Override
4621     public double getEstimatedP() {
4622         return estimatedP;
4623     }
4624 
4625     /**
4626      * Gets estimated measure of quality of estimated fit as a value between 0.0 and 1.0. The larger the quality value
4627      * is, the better the fit that has been estimated.
4628      *
4629      * @return estimated measure of quality of estimated fit.
4630      */
4631     @Override
4632     public double getEstimatedQ() {
4633         return estimatedQ;
4634     }
4635 
4636     /**
4637      * Gets variance of estimated x coordinate of gyroscope bias expressed in (rad^2/s^2).
4638      *
4639      * @return variance of estimated x coordinate of gyroscope bias or null if not available.
4640      */
4641     public Double getEstimatedBiasXVariance() {
4642         return estimatedCovariance != null ? estimatedCovariance.getElementAt(0, 0) : null;
4643     }
4644 
4645     /**
4646      * Gets standard deviation of estimated x coordinate of gyroscope bias expressed in
4647      * radians per second (rad/s).
4648      *
4649      * @return standard deviation of estimated x coordinate of gyroscope bias or null if not
4650      * available.
4651      */
4652     public Double getEstimatedBiasXStandardDeviation() {
4653         final var variance = getEstimatedBiasXVariance();
4654         return variance != null ? Math.sqrt(variance) : null;
4655     }
4656 
4657     /**
4658      * Gets standard deviation of estimated x coordinate of gyroscope bias.
4659      *
4660      * @return standard deviation of estimated x coordinate of gyroscope bias or null if not
4661      * available.
4662      */
4663     public AngularSpeed getEstimatedBiasXStandardDeviationAsAngularSpeed() {
4664         return estimatedCovariance != null
4665                 ? new AngularSpeed(getEstimatedBiasXStandardDeviation(), AngularSpeedUnit.RADIANS_PER_SECOND)
4666                 : null;
4667     }
4668 
4669     /**
4670      * Gets standard deviation of estimated x coordinate of gyroscope bias.
4671      *
4672      * @param result instance where result will be stored.
4673      * @return true if standard deviation of estimated x coordinate of gyroscope bias is available,
4674      * false otherwise.
4675      */
4676     public boolean getEstimatedBiasXStandardDeviationAsAngularSpeed(final AngularSpeed result) {
4677         if (estimatedCovariance != null) {
4678             result.setValue(getEstimatedBiasXStandardDeviation());
4679             result.setUnit(AngularSpeedUnit.RADIANS_PER_SECOND);
4680             return true;
4681         } else {
4682             return false;
4683         }
4684     }
4685 
4686     /**
4687      * Gets variance of estimated y coordinate of gyroscope bias expressed in (rad^2/s^2).
4688      *
4689      * @return variance of estimated y coordinate of gyroscope bias or null if not available.
4690      */
4691     public Double getEstimatedBiasYVariance() {
4692         return estimatedCovariance != null ? estimatedCovariance.getElementAt(1, 1) : null;
4693     }
4694 
4695     /**
4696      * Gets standard deviation of estimated y coordinate of gyroscope bias expressed in
4697      * radians per second (rad/s).
4698      *
4699      * @return standard deviation of estimated y coordinate of gyroscope bias or null if not
4700      * available.
4701      */
4702     public Double getEstimatedBiasYStandardDeviation() {
4703         final var variance = getEstimatedBiasYVariance();
4704         return variance != null ? Math.sqrt(variance) : null;
4705     }
4706 
4707     /**
4708      * Gets standard deviation of estimated y coordinate of gyroscope bias.
4709      *
4710      * @return standard deviation of estimated y coordinate of gyroscope bias or null if not
4711      * available.
4712      */
4713     public AngularSpeed getEstimatedBiasYStandardDeviationAsAngularSpeed() {
4714         return estimatedCovariance != null
4715                 ? new AngularSpeed(getEstimatedBiasYStandardDeviation(), AngularSpeedUnit.RADIANS_PER_SECOND)
4716                 : null;
4717     }
4718 
4719     /**
4720      * Gets standard deviation of estimated y coordinate of gyroscope bias.
4721      *
4722      * @param result instance where result will be stored.
4723      * @return true if standard deviation of estimated y coordinate of gyroscope bias is available,
4724      * false otherwise.
4725      */
4726     public boolean getEstimatedBiasYStandardDeviationAsAngularSpeed(final AngularSpeed result) {
4727         if (estimatedCovariance != null) {
4728             result.setValue(getEstimatedBiasYStandardDeviation());
4729             result.setUnit(AngularSpeedUnit.RADIANS_PER_SECOND);
4730             return true;
4731         } else {
4732             return false;
4733         }
4734     }
4735 
4736     /**
4737      * Gets variance of estimated z coordinate of gyroscope bias expressed in (rad^2/s^2).
4738      *
4739      * @return variance of estimated z coordinate of gyroscope bias or null if not available.
4740      */
4741     public Double getEstimatedBiasZVariance() {
4742         return estimatedCovariance != null ? estimatedCovariance.getElementAt(2, 2) : null;
4743     }
4744 
4745     /**
4746      * Gets standard deviation of estimated z coordinate of gyroscope bias expressed in
4747      * radians per second (rad/s).
4748      *
4749      * @return standard deviation of estimated z coordinate of gyroscope bias or null if not
4750      * available.
4751      */
4752     public Double getEstimatedBiasZStandardDeviation() {
4753         final var variance = getEstimatedBiasZVariance();
4754         return variance != null ? Math.sqrt(variance) : null;
4755     }
4756 
4757     /**
4758      * Gets standard deviation of estimated z coordinate of gyroscope bias.
4759      *
4760      * @return standard deviation of estimated z coordinate of gyroscope bias or null if not
4761      * available.
4762      */
4763     public AngularSpeed getEstimatedBiasZStandardDeviationAsAngularSpeed() {
4764         return estimatedCovariance != null
4765                 ? new AngularSpeed(getEstimatedBiasZStandardDeviation(), AngularSpeedUnit.RADIANS_PER_SECOND)
4766                 : null;
4767     }
4768 
4769     /**
4770      * Gets standard deviation of estimated z coordinate of gyroscope bias.
4771      *
4772      * @param result instance where result will be stored.
4773      * @return true if standard deviation of estimated z coordinate of gyroscope bias is available,
4774      * false otherwise.
4775      */
4776     public boolean getEstimatedBiasZStandardDeviationAsAngularSpeed(final AngularSpeed result) {
4777         if (estimatedCovariance != null) {
4778             result.setValue(getEstimatedBiasZStandardDeviation());
4779             result.setUnit(AngularSpeedUnit.RADIANS_PER_SECOND);
4780             return true;
4781         } else {
4782             return false;
4783         }
4784     }
4785 
4786     /**
4787      * Gets standard deviation of estimated gyroscope bias coordinates.
4788      *
4789      * @return standard deviation of estimated gyroscope bias coordinates.
4790      */
4791     public AngularSpeedTriad getEstimatedBiasStandardDeviation() {
4792         return estimatedCovariance != null
4793                 ? new AngularSpeedTriad(AngularSpeedUnit.RADIANS_PER_SECOND,
4794                 getEstimatedBiasXStandardDeviation(),
4795                 getEstimatedBiasYStandardDeviation(),
4796                 getEstimatedBiasZStandardDeviation())
4797                 : null;
4798     }
4799 
4800     /**
4801      * Gets standard deviation of estimated gyroscope bias coordinates.
4802      *
4803      * @param result instance where result will be stored.
4804      * @return true if standard deviation of gyroscope bias was available, false
4805      * otherwise.
4806      */
4807     public boolean getEstimatedBiasStandardDeviation(final AngularSpeedTriad result) {
4808         if (estimatedCovariance != null) {
4809             result.setValueCoordinatesAndUnit(
4810                     getEstimatedBiasXStandardDeviation(),
4811                     getEstimatedBiasYStandardDeviation(),
4812                     getEstimatedBiasZStandardDeviation(),
4813                     AngularSpeedUnit.RADIANS_PER_SECOND);
4814             return true;
4815         } else {
4816             return false;
4817         }
4818     }
4819 
4820     /**
4821      * Gets average of estimated standard deviation of gyroscope bias coordinates expressed
4822      * in radians per second (rad/s).
4823      *
4824      * @return average of estimated standard deviation of gyroscope bias coordinates or null
4825      * if not available.
4826      */
4827     public Double getEstimatedBiasStandardDeviationAverage() {
4828         return estimatedCovariance != null
4829                 ? (getEstimatedBiasXStandardDeviation() + getEstimatedBiasYStandardDeviation()
4830                 + getEstimatedBiasZStandardDeviation()) / 3.0
4831                 : null;
4832     }
4833 
4834     /**
4835      * Gets average of estimated standard deviation of gyroscope bias coordinates.
4836      *
4837      * @return average of estimated standard deviation of gyroscope bias coordinates or null.
4838      */
4839     public AngularSpeed getEstimatedBiasStandardDeviationAverageAsAngularSpeed() {
4840         return estimatedCovariance != null
4841                 ? new AngularSpeed(getEstimatedBiasStandardDeviationAverage(), AngularSpeedUnit.RADIANS_PER_SECOND)
4842                 : null;
4843     }
4844 
4845     /**
4846      * Gets average of estimated standard deviation of gyroscope bias coordinates.
4847      *
4848      * @param result instance where result will be stored.
4849      * @return true if average of estimated standard deviation of gyroscope bias is available,
4850      * false otherwise.
4851      */
4852     public boolean getEstimatedBiasStandardDeviationAverageAsAngularSpeed(final AngularSpeed result) {
4853         if (estimatedCovariance != null) {
4854             result.setValue(getEstimatedBiasStandardDeviationAverage());
4855             result.setUnit(AngularSpeedUnit.RADIANS_PER_SECOND);
4856             return true;
4857         } else {
4858             return false;
4859         }
4860     }
4861 
4862     /**
4863      * Gets norm of estimated standard deviation of gyroscope bias expressed in
4864      * radians per second (rad/s).
4865      * This can be used as the initial gyroscope bias uncertainty for
4866      * {@link INSLooselyCoupledKalmanInitializerConfig} or {@link INSTightlyCoupledKalmanInitializerConfig}.
4867      *
4868      * @return norm of estimated standard deviation of gyroscope bias or null
4869      * if not available.
4870      */
4871     @Override
4872     public Double getEstimatedBiasStandardDeviationNorm() {
4873         return estimatedCovariance != null
4874                 ? Math.sqrt(getEstimatedBiasXVariance() + getEstimatedBiasYVariance() + getEstimatedBiasZVariance())
4875                 : null;
4876     }
4877 
4878     /**
4879      * Gets norm of estimated standard deviation of gyroscope bias.
4880      * This can be used as the initial gyroscope bias uncertainty for
4881      * {@link INSLooselyCoupledKalmanInitializerConfig} or {@link INSTightlyCoupledKalmanInitializerConfig}.
4882      *
4883      * @return norm of estimated standard deviation of gyroscope bias or null
4884      * if not available.
4885      */
4886     public AngularSpeed getEstimatedBiasStandardDeviationNormAsAngularSpeed() {
4887         return estimatedCovariance != null
4888                 ? new AngularSpeed(getEstimatedBiasStandardDeviationNorm(), AngularSpeedUnit.RADIANS_PER_SECOND)
4889                 : null;
4890     }
4891 
4892     /**
4893      * Gets norm of estimated standard deviation of gyroscope bias coordinates.
4894      * This can be used as the initial gyroscope bias uncertainty for
4895      * {@link INSLooselyCoupledKalmanInitializerConfig} or {@link INSTightlyCoupledKalmanInitializerConfig}.
4896      *
4897      * @param result instance where result will be stored.
4898      * @return true if norm of estimated standard deviation of gyroscope bias is
4899      * available, false otherwise.
4900      */
4901     public boolean getEstimatedBiasStandardDeviationNormAsAngularSpeed(final AngularSpeed result) {
4902         if (estimatedCovariance != null) {
4903             result.setValue(getEstimatedBiasStandardDeviationNorm());
4904             result.setUnit(AngularSpeedUnit.RADIANS_PER_SECOND);
4905             return true;
4906         } else {
4907             return false;
4908         }
4909     }
4910 
4911     /**
4912      * Internal method to perform calibration when common z-axis is assumed
4913      * for both the accelerometer and gyroscope and when G-dependent cross
4914      * biases are being estimated.
4915      *
4916      * @throws AlgebraException                              if there are numerical errors.
4917      * @throws FittingException                              if no convergence to solution is found.
4918      * @throws com.irurueta.numerical.NotReadyException      if fitter is not ready.
4919      * @throws InvalidSourceAndDestinationFrameTypeException never happens
4920      */
4921     private void calibrateCommonAxisAndGDependentCrossBiases() throws AlgebraException, FittingException,
4922             com.irurueta.numerical.NotReadyException, InvalidSourceAndDestinationFrameTypeException {
4923 
4924         // The gyroscope model is
4925         // Ωmeas = bg + (I + Mg) * Ωtrue + Gg * ftrue + w
4926 
4927         // Ideally a least squares solution tries to minimize noise component, so:
4928         // Ωmeas = bg + (I + Mg) * Ωtrue + Gg * ftrue
4929 
4930         // For convergence purposes of the Levenberg-Marquardt algorithm, we
4931         // take common factor M = I + Mg
4932 
4933         // and the gyroscope model can be better expressed as:
4934 
4935         // Ωmeas = M*(Ωtrue + b + G * ftrue)
4936 
4937         // where:
4938         // bg = M*b --> b = M^-1*bg
4939         // Gg = M*G --> G = M^-1*Gg
4940 
4941         // We know that the norm of the true angular rate when the device is in a fixed
4942         // and unknown position and orientation is equal to the Earth rotation rate.
4943         // ||Ωtrue|| = 7.292115E-5 rad/s
4944 
4945         // Hence
4946         // Ωmeas - M*b - M*G*ftrue = M*Ωtrue
4947         // M^-1 * (Ωmeas - M*b - M*G*ftrue) = Ωtrue
4948 
4949         // ||Ωtrue||^2 = (M^-1 * (Ωmeas - M*b - M*G*ftrue))^T*(M^-1 * (Ωmeas - M*b - M*G*ftrue))
4950         // ||Ωtrue||^2 = (Ωmeas - M*b - M*G*ftrue)^T * (M^-1)^T * M^-1 * (Ωmeas - M*b - M*G*ftrue)
4951         // ||Ωtrue||^2 = (Ωmeas - M*b - M*G*ftrue)^T * ||M^-1||^2 * (Ωmeas - M*b - M*G*ftrue)
4952         // ||Ωtrue||^2 = ||Ωmeas - M*b - M*G*ftrue||^2 * ||M^-1||^2
4953 
4954         // Where:
4955 
4956         // b = [bx]
4957         //     [by]
4958         //     [bz]
4959 
4960         // M = [m11 	m12 	m13]
4961         //     [0 		m22 	m23]
4962         //     [0 	 	0 		m33]
4963 
4964         // G = [g11 	g12 	g13]
4965         //     [g21 	g22 	g23]
4966         //     [g31 	g32 	g33]
4967 
4968         // ftrue = [ftruex]
4969         //         [ftruey]
4970         //         [fturez]
4971 
4972         final var gradientEstimator = new GradientEstimator(this::evaluateCommonAxisWithGDependentCrossBiases);
4973 
4974         final var initialM = Matrix.identity(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
4975         initialM.add(getInitialMg());
4976 
4977         // Force initial M to be upper diagonal
4978         initialM.setElementAt(1, 0, 0.0);
4979         initialM.setElementAt(2, 0, 0.0);
4980         initialM.setElementAt(2, 1, 0.0);
4981 
4982         final var invInitM = Utils.inverse(initialM);
4983         final var initBg = getInitialBiasAsMatrix();
4984         final var initB = invInitM.multiplyAndReturnNew(initBg);
4985         final var initGg = getInitialGg();
4986         final var initG = invInitM.multiplyAndReturnNew(initGg);
4987 
4988         fitter.setFunctionEvaluator(new LevenbergMarquardtMultiDimensionFunctionEvaluator() {
4989             @Override
4990             public int getNumberOfDimensions() {
4991                 // Input points are measured angular rate coordinates +
4992                 // measured specific force coordinates
4993                 return 2 * BodyKinematics.COMPONENTS;
4994             }
4995 
4996             @Override
4997             public double[] createInitialParametersArray() {
4998                 final var initial = new double[COMMON_Z_AXIS_UNKNOWNS_AND_CROSS_BIASES];
4999 
5000                 // biases b
5001                 for (var i = 0; i < BodyKinematics.COMPONENTS; i++) {
5002                     initial[i] = initB.getElementAtIndex(i);
5003                 }
5004 
5005                 // upper diagonal cross coupling errors M
5006                 var k = BodyKinematics.COMPONENTS;
5007                 for (var j = 0; j < BodyKinematics.COMPONENTS; j++) {
5008                     for (var i = 0; i < BodyKinematics.COMPONENTS; i++) {
5009                         if (i <= j) {
5010                             initial[k] = initialM.getElementAt(i, j);
5011                             k++;
5012                         }
5013                     }
5014                 }
5015 
5016                 // g-dependent cross biases G
5017                 final var num = BodyKinematics.COMPONENTS * BodyKinematics.COMPONENTS;
5018                 for (int i = 0, j = k; i < num; i++, j++) {
5019                     initial[j] = initG.getElementAtIndex(i);
5020                 }
5021 
5022                 return initial;
5023             }
5024 
5025             @Override
5026             public double evaluate(
5027                     final int i, final double[] point, final double[] params, final double[] derivatives)
5028                     throws EvaluationException {
5029 
5030                 measAngularRateX = point[0];
5031                 measAngularRateY = point[1];
5032                 measAngularRateZ = point[2];
5033 
5034                 fmeasX = point[3];
5035                 fmeasY = point[4];
5036                 fmeasZ = point[5];
5037 
5038                 gradientEstimator.gradient(params, derivatives);
5039 
5040                 return evaluateCommonAxisWithGDependentCrossBiases(params);
5041             }
5042         });
5043 
5044         setInputDataWithGDependentCrossBiases();
5045 
5046         fitter.fit();
5047 
5048         final var result = fitter.getA();
5049 
5050         final var bx = result[0];
5051         final var by = result[1];
5052         final var bz = result[2];
5053 
5054         final var m11 = result[3];
5055 
5056         final var m12 = result[4];
5057         final var m22 = result[5];
5058 
5059         final var m13 = result[6];
5060         final var m23 = result[7];
5061         final var m33 = result[8];
5062 
5063         final var g11 = result[9];
5064         final var g21 = result[10];
5065         final var g31 = result[11];
5066 
5067         final var g12 = result[12];
5068         final var g22 = result[13];
5069         final var g32 = result[14];
5070 
5071         final var g13 = result[15];
5072         final var g23 = result[16];
5073         final var g33 = result[17];
5074 
5075         final var mb = new Matrix(BodyKinematics.COMPONENTS, 1);
5076         mb.setElementAtIndex(0, bx);
5077         mb.setElementAtIndex(1, by);
5078         mb.setElementAtIndex(2, bz);
5079 
5080         final var mm = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
5081         mm.setElementAtIndex(0, m11);
5082         mm.setElementAtIndex(1, 0.0);
5083         mm.setElementAtIndex(2, 0.0);
5084 
5085         mm.setElementAtIndex(3, m12);
5086         mm.setElementAtIndex(4, m22);
5087         mm.setElementAtIndex(5, 0.0);
5088 
5089         mm.setElementAtIndex(6, m13);
5090         mm.setElementAtIndex(7, m23);
5091         mm.setElementAtIndex(8, m33);
5092 
5093         final var mg = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
5094         mg.setElementAtIndex(0, g11);
5095         mg.setElementAtIndex(1, g21);
5096         mg.setElementAtIndex(2, g31);
5097 
5098         mg.setElementAtIndex(3, g12);
5099         mg.setElementAtIndex(4, g22);
5100         mg.setElementAtIndex(5, g32);
5101 
5102         mg.setElementAtIndex(6, g13);
5103         mg.setElementAtIndex(7, g23);
5104         mg.setElementAtIndex(8, g33);
5105 
5106         setResult(mm, mb, mg);
5107 
5108         // at this point covariance is expressed in terms of b, M and G, and must
5109         // be expressed in terms of bg, Mg and Gg.
5110         // We know that:
5111         // bg = M * b
5112         // Mg = M - I
5113         // Gg = M * G
5114 
5115         // b = [bx]
5116         //     [by]
5117         //     [bz]
5118 
5119         // M = [m11  m12  m13]
5120         //     [0    m22  m23]
5121         //     [0    0    m33]
5122 
5123         // G = [g11  g12  g13]
5124         //     [g21  g22  g23]
5125         //     [g31  g32  g33]
5126 
5127         // bg = [m11  m12  m13][bx] = [m11 * bx + m12 * by + m13 * bz] = [bgx]
5128         //      [0    m22  m23][by]   [           m22 * by + m23 * bz]   [bgy]
5129         //      [0    0    m33][bz]   [                      m33 * bz]   [bgz]
5130 
5131         // Mg = [sx  mxy  mxz] = [m11 - 1    m12         m13     ]
5132         //      [myx sy	  myz]   [0          m22 - 1     m23     ]
5133         //      [mzx mzy  sz ]   [0          0           m33 - 1 ]
5134 
5135         // Gg = [gg11  gg12  gg13] = [m11  m12  m13][g11  g12  g13]
5136         //      [gg21  gg22  gg23]   [0    m22  m23][g21  g22  g23]
5137         //      [gg31  gg32  gg33]   [0    0    m33][g31  g32  g33]
5138 
5139         // Defining the linear application:
5140         // F(b, M, G) = F(bx, by, bz, m11, m12, m22, m13, m23, m33, g11, g21, g31, g12, g22, g32, g13, g23, g33)
5141         // as:
5142         // [bgx] =  [m11 * bx + m12 * by + m13 * bz]
5143         // [bgy]    [           m22 * by + m23 * bz]
5144         // [bgz]    [                      m33 * bz]
5145         // [sx]     [m11 - 1]
5146         // [sy]     [m22 - 1]
5147         // [sz]     [m33 - 1]
5148         // [mxy]    [m12]
5149         // [mxz]    [m13]
5150         // [myx]    [0]
5151         // [myz]    [m23]
5152         // [mzx]    [0]
5153         // [mzy]    [0]
5154         // [gg11]   [m11 * g11 + m12 * g21 + m13 * g31]
5155         // [gg21]   [            m22 * g21 + m23 * g31]
5156         // [gg31]   [                        m33 * g31]
5157         // [gg12]   [m11 * g12 + m12 * g22 + m13 * g32]
5158         // [gg22]   [            m22 * g22 + m23 * g32]
5159         // [gg32]   [                        m33 * g32]
5160         // [gg13]   [m11 * g13 + m12 * g23 + m13 * g33]
5161         // [gg23]   [            m22 * g23 + m23 * g33]
5162         // [gg33]   [                        m33 * g33]
5163 
5164         // Then the Jacobian of F(b, M, G) is:
5165         // J = [m11  m12  m13  bx   by   0    bz   0    0    0    0    0    0    0    0    0    0    0  ]
5166         //     [0    m22  m23  0    0    by   0    bz   0    0    0    0    0    0    0    0    0    0  ]
5167         //     [0    0    m33  0    0    0    0    0    bz   0    0    0    0    0    0    0    0    0  ]
5168         //     [0    0    0    1    0    0    0    0    0    0    0    0    0    0    0    0    0    0  ]
5169         //     [0    0    0    0    0    1    0    0    0    0    0    0    0    0    0    0    0    0  ]
5170         //     [0    0    0    0    0    0    0    0    1    0    0    0    0    0    0    0    0    0  ]
5171         //     [0    0    0    0    1    0    0    0    0    0    0    0    0    0    0    0    0    0  ]
5172         //     [0    0    0    0    0    0    1    0    0    0    0    0    0    0    0    0    0    0  ]
5173         //     [0    0    0    0    0    0    0    0    0    0    0    0    0    0    0    0    0    0  ]
5174         //     [0    0    0    0    0    0    0    1    0    0    0    0    0    0    0    0    0    0  ]
5175         //     [0    0    0    0    0    0    0    0    0    0    0    0    0    0    0    0    0    0  ]
5176         //     [0    0    0    0    0    0    0    0    0    0    0    0    0    0    0    0    0    0  ]
5177         //     [0    0    0    g11  g21  0    g31  0    0    m11  m12  m13  0    0    0    0    0    0  ]
5178         //     [0    0    0    0    0    g21  0    g31  0    0    m22  m23  0    0    0    0    0    0  ]
5179         //     [0    0    0    0    0    0    0    0    g31  0    0    m33  0    0    0    0    0    0  ]
5180         //     [0    0    0    g12  g22  0    g32  0    0    0    0    0    m11  m12  m13  0    0    0  ]
5181         //     [0    0    0    0    0    g22  0    g32  0    0    0    0    0    m22  m23  0    0    0  ]
5182         //     [0    0    0    0    0    0    0    0    g32  0    0    0    0    0    m33  0    0    0  ]
5183         //     [0    0    0    g13  g23  0    g33  0    0    0    0    0    0    0    0    m11  m12  m13]
5184         //     [0    0    0    0    0    g23  0    g33  0    0    0    0    0    0    0    0    m22  m23]
5185         //     [0    0    0    0    0    0    0    0    g33  0    0    0    0    0    0    0    0    m33]
5186 
5187         // We know that the propagated covariance is J * Cov * J', hence:
5188         final var jacobian = new Matrix(GENERAL_UNKNOWNS_AND_CROSS_BIASES, COMMON_Z_AXIS_UNKNOWNS_AND_CROSS_BIASES);
5189 
5190         jacobian.setElementAt(0, 0, m11);
5191         jacobian.setElementAt(0, 1, m12);
5192         jacobian.setElementAt(0, 2, m13);
5193         jacobian.setElementAt(0, 3, bx);
5194         jacobian.setElementAt(0, 4, by);
5195         jacobian.setElementAt(0, 6, bz);
5196 
5197         jacobian.setElementAt(1, 1, m22);
5198         jacobian.setElementAt(1, 2, m23);
5199         jacobian.setElementAt(1, 5, by);
5200         jacobian.setElementAt(1, 7, bz);
5201 
5202         jacobian.setElementAt(2, 2, m33);
5203         jacobian.setElementAt(2, 8, bz);
5204 
5205         jacobian.setElementAt(3, 3, 1.0);
5206         jacobian.setElementAt(4, 5, 1.0);
5207         jacobian.setElementAt(5, 8, 1.0);
5208 
5209         jacobian.setElementAt(6, 4, 1.0);
5210         jacobian.setElementAt(7, 6, 1.0);
5211 
5212         jacobian.setElementAt(9, 7, 1.0);
5213 
5214         jacobian.setElementAt(12, 3, g11);
5215         jacobian.setElementAt(12, 4, g21);
5216         jacobian.setElementAt(12, 6, g31);
5217         jacobian.setElementAt(12, 9, m11);
5218         jacobian.setElementAt(12, 10, m12);
5219         jacobian.setElementAt(12, 11, m13);
5220 
5221         jacobian.setElementAt(13, 5, g21);
5222         jacobian.setElementAt(13, 7, g31);
5223         jacobian.setElementAt(13, 10, m22);
5224         jacobian.setElementAt(13, 11, m23);
5225 
5226         jacobian.setElementAt(14, 8, g31);
5227         jacobian.setElementAt(14, 11, m33);
5228 
5229         jacobian.setElementAt(15, 3, g12);
5230         jacobian.setElementAt(15, 4, g22);
5231         jacobian.setElementAt(15, 6, g32);
5232         jacobian.setElementAt(15, 12, m11);
5233         jacobian.setElementAt(15, 13, m12);
5234         jacobian.setElementAt(15, 14, m13);
5235 
5236         jacobian.setElementAt(16, 5, g22);
5237         jacobian.setElementAt(16, 7, g32);
5238         jacobian.setElementAt(16, 13, m22);
5239         jacobian.setElementAt(16, 14, m23);
5240 
5241         jacobian.setElementAt(17, 8, g32);
5242         jacobian.setElementAt(17, 14, m33);
5243 
5244         jacobian.setElementAt(18, 3, g13);
5245         jacobian.setElementAt(18, 4, g23);
5246         jacobian.setElementAt(18, 6, g33);
5247         jacobian.setElementAt(18, 15, m11);
5248         jacobian.setElementAt(18, 16, m12);
5249         jacobian.setElementAt(18, 17, m13);
5250 
5251         jacobian.setElementAt(19, 5, g23);
5252         jacobian.setElementAt(19, 7, g33);
5253         jacobian.setElementAt(19, 16, m22);
5254         jacobian.setElementAt(19, 17, m23);
5255 
5256         jacobian.setElementAt(20, 8, g33);
5257         jacobian.setElementAt(20, 17, m33);
5258 
5259         final var jacobianTrans = jacobian.transposeAndReturnNew();
5260         jacobian.multiply(estimatedCovariance);
5261         jacobian.multiply(jacobianTrans);
5262         estimatedCovariance = jacobian;
5263     }
5264 
5265     /**
5266      * Internal method to perform general calibration when G-dependent cross
5267      * biases are being estimated.
5268      *
5269      * @throws AlgebraException                              if there are numerical errors.
5270      * @throws FittingException                              if no convergence to solution is found.
5271      * @throws com.irurueta.numerical.NotReadyException      if fitter is not ready.
5272      * @throws InvalidSourceAndDestinationFrameTypeException never happens
5273      */
5274     private void calibrateGeneralAndGDependentCrossBiases() throws AlgebraException, FittingException,
5275             com.irurueta.numerical.NotReadyException, InvalidSourceAndDestinationFrameTypeException {
5276 
5277         // The gyroscope model is
5278         // Ωmeas = bg + (I + Mg) * Ωtrue + Gg * ftrue + w
5279 
5280         // Ideally a least squares solution tries to minimize noise component, so:
5281         // Ωmeas = bg + (I + Mg) * Ωtrue + Gg * ftrue
5282 
5283         // For convergence purposes of the Levenberg-Marquardt algorithm, we
5284         // take common factor M = I + Mg
5285 
5286         // and the gyroscope model can be better expressed as:
5287 
5288         // Ωmeas = M*(Ωtrue + b + G * ftrue)
5289 
5290         // where:
5291         // bg = M*b --> b = M^-1*bg
5292         // Gg = M*G --> G = M^-1*Gg
5293 
5294         // We know that the norm of the true angular rate when the device is in a pixed
5295         // and unknown position and orientation is equal to the Earth rotation rate.
5296         // ||Ωtrue|| = 7.292115E-5 rad/s
5297 
5298         // Hence
5299         // Ωmeas - M*b - M*G*ftrue = M*Ωtrue
5300         // M^-1 * (Ωmeas - M*b - M*G*ftrue) = Ωtrue
5301 
5302         // ||Ωtrue||^2 = (M^-1 * (Ωmeas - M*b - M*G*ftrue))^T*(M^-1 * (Ωmeas - M*b - M*G*ftrue))
5303         // ||Ωtrue||^2 = (Ωmeas - M*b - M*G*ftrue)^T * (M^-1)^T * M^-1 * (Ωmeas - M*b - M*G*ftrue)
5304         // ||Ωtrue||^2 = (Ωmeas - M*b - M*G*ftrue)^T * ||M^-1||^2 * (Ωmeas - M*b - M*G*ftrue)
5305         // ||Ωtrue||^2 = ||Ωmeas - M*b - M*G*ftrue||^2 * ||M^-1||^2
5306 
5307         // Where:
5308 
5309         // b = [bx]
5310         //     [by]
5311         //     [bz]
5312 
5313         // M = [m11 	m12 	m13]
5314         //     [m21 	m22 	m23]
5315         //     [m31 	m32 	m33]
5316 
5317         // G = [g11 	g12 	g13]
5318         //     [g21 	g22 	g23]
5319         //     [g31 	g32 	g33]
5320 
5321         // ftrue = [ftruex]
5322         //         [ftruey]
5323         //         [fturez]
5324 
5325         final var gradientEstimator = new GradientEstimator(this::evaluateGeneralWithGDependentCrossBiases);
5326 
5327         final var initialM = Matrix.identity(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
5328         initialM.add(getInitialMg());
5329 
5330         final var invInitM = Utils.inverse(initialM);
5331         final var initBg = getInitialBiasAsMatrix();
5332         final var initB = invInitM.multiplyAndReturnNew(initBg);
5333         final var initGg = getInitialGg();
5334         final var initG = invInitM.multiplyAndReturnNew(initGg);
5335 
5336         fitter.setFunctionEvaluator(new LevenbergMarquardtMultiDimensionFunctionEvaluator() {
5337             @Override
5338             public int getNumberOfDimensions() {
5339                 // Input points are measured angular rate coordinates +
5340                 // measured specific force coordinates
5341                 return 2 * BodyKinematics.COMPONENTS;
5342             }
5343 
5344             @Override
5345             public double[] createInitialParametersArray() {
5346                 final var initial = new double[GENERAL_UNKNOWNS_AND_CROSS_BIASES];
5347 
5348                 // biases b
5349                 for (var i = 0; i < BodyKinematics.COMPONENTS; i++) {
5350                     initial[i] = initB.getElementAtIndex(i);
5351                 }
5352 
5353                 // cross coupling errors M
5354                 final var num = BodyKinematics.COMPONENTS * BodyKinematics.COMPONENTS;
5355                 for (int i = 0, j = BodyKinematics.COMPONENTS; i < num; i++, j++) {
5356                     initial[j] = initialM.getElementAtIndex(i);
5357                 }
5358 
5359                 // g-dependent cross biases G
5360                 for (int i = 0, j = BodyKinematics.COMPONENTS + num; i < num; i++, j++) {
5361                     initial[j] = initG.getElementAtIndex(i);
5362                 }
5363 
5364                 return initial;
5365             }
5366 
5367             @Override
5368             public double evaluate(
5369                     final int i, final double[] point, final double[] params, final double[] derivatives)
5370                     throws EvaluationException {
5371 
5372                 measAngularRateX = point[0];
5373                 measAngularRateY = point[1];
5374                 measAngularRateZ = point[2];
5375 
5376                 fmeasX = point[3];
5377                 fmeasY = point[4];
5378                 fmeasZ = point[5];
5379 
5380                 gradientEstimator.gradient(params, derivatives);
5381 
5382                 return evaluateGeneralWithGDependentCrossBiases(params);
5383             }
5384         });
5385 
5386         setInputDataWithGDependentCrossBiases();
5387 
5388         fitter.fit();
5389 
5390         final var result = fitter.getA();
5391 
5392         final var bx = result[0];
5393         final var by = result[1];
5394         final var bz = result[2];
5395 
5396         final var m11 = result[3];
5397         final var m21 = result[4];
5398         final var m31 = result[5];
5399 
5400         final var m12 = result[6];
5401         final var m22 = result[7];
5402         final var m32 = result[8];
5403 
5404         final var m13 = result[9];
5405         final var m23 = result[10];
5406         final var m33 = result[11];
5407 
5408         final var g11 = result[12];
5409         final var g21 = result[13];
5410         final var g31 = result[14];
5411 
5412         final var g12 = result[15];
5413         final var g22 = result[16];
5414         final var g32 = result[17];
5415 
5416         final var g13 = result[18];
5417         final var g23 = result[19];
5418         final var g33 = result[20];
5419 
5420         final var mb = new Matrix(BodyKinematics.COMPONENTS, 1);
5421         mb.setElementAtIndex(0, bx);
5422         mb.setElementAtIndex(1, by);
5423         mb.setElementAtIndex(2, bz);
5424 
5425         final var mm = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
5426         mm.setElementAtIndex(0, m11);
5427         mm.setElementAtIndex(1, m21);
5428         mm.setElementAtIndex(2, m31);
5429 
5430         mm.setElementAtIndex(3, m12);
5431         mm.setElementAtIndex(4, m22);
5432         mm.setElementAtIndex(5, m32);
5433 
5434         mm.setElementAtIndex(6, m13);
5435         mm.setElementAtIndex(7, m23);
5436         mm.setElementAtIndex(8, m33);
5437 
5438         final var mg = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
5439         mg.setElementAtIndex(0, g11);
5440         mg.setElementAtIndex(1, g21);
5441         mg.setElementAtIndex(2, g31);
5442 
5443         mg.setElementAtIndex(3, g12);
5444         mg.setElementAtIndex(4, g22);
5445         mg.setElementAtIndex(5, g32);
5446 
5447         mg.setElementAtIndex(6, g13);
5448         mg.setElementAtIndex(7, g23);
5449         mg.setElementAtIndex(8, g33);
5450 
5451         setResult(mm, mb, mg);
5452 
5453         // at this point covariance is expressed in terms of b, M and G, and must
5454         // be expressed in terms of bg, Mg and Gg.
5455         // We know that:
5456         // bg = M * b
5457         // Mg = M - I
5458         // Gg = M * G
5459 
5460         // b = [bx]
5461         //     [by]
5462         //     [bz]
5463 
5464         // M = [m11  m12  m13]
5465         //     [m21  m22  m23]
5466         //     [m31  m32  m33]
5467 
5468         // G = [g11  g12  g13]
5469         //     [g21  g22  g23]
5470         //     [g31  g32  g33]
5471 
5472         // bg = [m11  m12  m13][bx] = [m11 * bx + m12 * by + m13 * bz] = [bgx]
5473         //      [m21  m22  m23][by]   [m21 * bx + m22 * by + m23 * bz]   [bgy]
5474         //      [m31  m32  m33][bz]   [m31 * bx + m32 * by + m33 * bz]   [bgz]
5475 
5476         // Mg = [sx  mxy  mxz] = [m11 - 1    m12         m13     ]
5477         //      [myx sy	  myz]   [m21        m22 - 1     m23     ]
5478         //      [mzx mzy  sz ]   [m31        m32         m33 - 1 ]
5479 
5480         // Gg = [gg11  gg12  gg13] = [m11  m12  m13][g11  g12  g13]
5481         //      [gg21  gg22  gg23]   [m21  m22  m23][g21  g22  g23]
5482         //      [gg31  gg32  gg33]   [m31  m32  m33][g31  g32  g33]
5483 
5484         // Defining the linear application:
5485         // F(b, M, G) = F(bx, by, bz, m11, m21, m31, m12, m22, m32, m13, m23, m33, g11, g21, g31, g12, g22, g32, g13, g23, g33)
5486         // as:
5487         // [bgx] =  [m11 * bx + m12 * by + m13 * bz]
5488         // [bgy]    [m21 * bx + m22 * by + m23 * bz]
5489         // [bgz]    [m31 * bx + m32 * by + m33 * bz]
5490         // [sx]     [m11 - 1]
5491         // [sy]     [m22 - 1]
5492         // [sz]     [m33 - 1]
5493         // [mxy]    [m12]
5494         // [mxz]    [m13]
5495         // [myx]    [m21]
5496         // [myz]    [m23]
5497         // [mzx]    [m31]
5498         // [mzy]    [m32]
5499         // [gg11]   [m11 * g11 + m12 * g21 + m13 * g31]
5500         // [gg21]   [m21 * g11 + m22 * g21 + m23 * g31]
5501         // [gg31]   [m31 * g11 + m32 * g21 + m33 * g31]
5502         // [gg12]   [m11 * g12 + m12 * g22 + m13 * g32]
5503         // [gg22]   [m21 * g12 + m22 * g22 + m23 * g32]
5504         // [gg32]   [m31 * g12 + m32 * g22 + m33 * g32]
5505         // [gg13]   [m11 * g13 + m12 * g23 + m13 * g33]
5506         // [gg23]   [m21 * g13 + m22 * g23 + m23 * g33]
5507         // [gg33]   [m31 * g13 + m32 * g23 + m33 * g33]
5508 
5509         // Then the Jacobian of F(b, M, G) is:
5510         // J = [m11  m12  m13  bx   0    0    by   0    0    bz   0    0    0    0    0    0    0    0    0    0    0  ]
5511         //     [m21  m22  m23  0    bx   0    0    by   0    0    bz   0    0    0    0    0    0    0    0    0    0  ]
5512         //     [m31  m32  m33  0    0    bx   0    0    by   0    0    bz   0    0    0    0    0    0    0    0    0  ]
5513         //     [0    0    0    1    0    0    0    0    0    0    0    0    0    0    0    0    0    0    0    0    0  ]
5514         //     [0    0    0    0    0    0    0    1    0    0    0    0    0    0    0    0    0    0    0    0    0  ]
5515         //     [0    0    0    0    0    0    0    0    0    0    0    1    0    0    0    0    0    0    0    0    0  ]
5516         //     [0    0    0    0    0    0    1    0    0    0    0    0    0    0    0    0    0    0    0    0    0  ]
5517         //     [0    0    0    0    0    0    0    0    0    1    0    0    0    0    0    0    0    0    0    0    0  ]
5518         //     [0    0    0    0    1    0    0    0    0    0    0    0    0    0    0    0    0    0    0    0    0  ]
5519         //     [0    0    0    0    0    0    0    0    0    0    1    0    0    0    0    0    0    0    0    0    0  ]
5520         //     [0    0    0    0    0    1    0    0    0    0    0    0    0    0    0    0    0    0    0    0    0  ]
5521         //     [0    0    0    0    0    0    0    0    1    0    0    0    0    0    0    0    0    0    0    0    0  ]
5522         //     [0    0    0    g11  0    0    g21  0    0    g31  0    0    m11  m12  m13  0    0    0    0    0    0  ]
5523         //     [0    0    0    0    g11  0    0    g21  0    0    g31  0    m21  m22  m23  0    0    0    0    0    0  ]
5524         //     [0    0    0    0    0    g11  0    0    g21  0    0    g31  m31  m32  m33  0    0    0    0    0    0  ]
5525         //     [0    0    0    g12  0    0    g22  0    0    g32  0    0    0    0    0    m11  m12  m13  0    0    0  ]
5526         //     [0    0    0    0    g12  0    0    g22  0    0    g32  0    0    0    0    m21  m22  m23  0    0    0  ]
5527         //     [0    0    0    0    0    g12  0    0    g22  0    0    g32  0    0    0    m31  m32  m33  0    0    0  ]
5528         //     [0    0    0    g13  0    0    g23  0    0    g33  0    0    0    0    0    0    0    0    m11  m12  m13]
5529         //     [0    0    0    0    g13  0    0    g23  0    0    g33  0    0    0    0    0    0    0    m21  m22  m23]
5530         //     [0    0    0    0    0    g13  0    0    g23  0    0    g33  0    0    0    0    0    0    m31  m32  m33]
5531 
5532         // We know that the propagated covariance is J * Cov * J', hence:
5533         final var jacobian = new Matrix(GENERAL_UNKNOWNS_AND_CROSS_BIASES, GENERAL_UNKNOWNS_AND_CROSS_BIASES);
5534 
5535         jacobian.setElementAt(0, 0, m11);
5536         jacobian.setElementAt(0, 1, m12);
5537         jacobian.setElementAt(0, 2, m13);
5538         jacobian.setElementAt(0, 3, bx);
5539         jacobian.setElementAt(0, 6, by);
5540         jacobian.setElementAt(0, 9, bz);
5541 
5542         jacobian.setElementAt(1, 0, m21);
5543         jacobian.setElementAt(1, 1, m22);
5544         jacobian.setElementAt(1, 2, m23);
5545         jacobian.setElementAt(1, 4, bx);
5546         jacobian.setElementAt(1, 7, by);
5547         jacobian.setElementAt(1, 10, bz);
5548 
5549         jacobian.setElementAt(2, 0, m31);
5550         jacobian.setElementAt(2, 1, m32);
5551         jacobian.setElementAt(2, 2, m33);
5552         jacobian.setElementAt(2, 5, bx);
5553         jacobian.setElementAt(2, 8, by);
5554         jacobian.setElementAt(2, 11, bz);
5555 
5556         jacobian.setElementAt(3, 3, 1.0);
5557         jacobian.setElementAt(4, 7, 1.0);
5558         jacobian.setElementAt(5, 11, 1.0);
5559 
5560         jacobian.setElementAt(6, 6, 1.0);
5561         jacobian.setElementAt(7, 9, 1.0);
5562         jacobian.setElementAt(8, 4, 1.0);
5563 
5564         jacobian.setElementAt(9, 10, 1.0);
5565         jacobian.setElementAt(10, 5, 1.0);
5566         jacobian.setElementAt(11, 8, 1.9);
5567 
5568         jacobian.setElementAt(12, 3, g11);
5569         jacobian.setElementAt(12, 6, g21);
5570         jacobian.setElementAt(12, 9, g31);
5571         jacobian.setElementAt(12, 12, m11);
5572         jacobian.setElementAt(12, 13, m12);
5573         jacobian.setElementAt(12, 14, m13);
5574 
5575         jacobian.setElementAt(13, 4, g11);
5576         jacobian.setElementAt(13, 7, g21);
5577         jacobian.setElementAt(13, 10, g31);
5578         jacobian.setElementAt(13, 12, m21);
5579         jacobian.setElementAt(13, 13, m22);
5580         jacobian.setElementAt(13, 14, m23);
5581 
5582         jacobian.setElementAt(14, 5, g11);
5583         jacobian.setElementAt(14, 8, g21);
5584         jacobian.setElementAt(14, 11, g31);
5585         jacobian.setElementAt(14, 12, m31);
5586         jacobian.setElementAt(14, 13, m32);
5587         jacobian.setElementAt(14, 14, m33);
5588 
5589         jacobian.setElementAt(15, 3, g12);
5590         jacobian.setElementAt(15, 6, g22);
5591         jacobian.setElementAt(15, 9, g32);
5592         jacobian.setElementAt(15, 15, m11);
5593         jacobian.setElementAt(15, 16, m12);
5594         jacobian.setElementAt(15, 17, m13);
5595 
5596         jacobian.setElementAt(16, 4, g12);
5597         jacobian.setElementAt(16, 7, g22);
5598         jacobian.setElementAt(16, 10, g32);
5599         jacobian.setElementAt(16, 15, m21);
5600         jacobian.setElementAt(16, 16, m22);
5601         jacobian.setElementAt(16, 17, m23);
5602 
5603         jacobian.setElementAt(17, 5, g12);
5604         jacobian.setElementAt(17, 8, g22);
5605         jacobian.setElementAt(17, 11, g32);
5606         jacobian.setElementAt(17, 15, m31);
5607         jacobian.setElementAt(17, 16, m32);
5608         jacobian.setElementAt(17, 17, m33);
5609 
5610         jacobian.setElementAt(18, 3, g13);
5611         jacobian.setElementAt(18, 6, g23);
5612         jacobian.setElementAt(18, 9, g33);
5613         jacobian.setElementAt(18, 18, m11);
5614         jacobian.setElementAt(18, 19, m12);
5615         jacobian.setElementAt(18, 20, m13);
5616 
5617         jacobian.setElementAt(19, 4, g13);
5618         jacobian.setElementAt(19, 7, g23);
5619         jacobian.setElementAt(19, 10, g33);
5620         jacobian.setElementAt(19, 18, m21);
5621         jacobian.setElementAt(19, 19, m22);
5622         jacobian.setElementAt(19, 20, m23);
5623 
5624         jacobian.setElementAt(20, 5, g13);
5625         jacobian.setElementAt(20, 8, g23);
5626         jacobian.setElementAt(20, 11, g33);
5627         jacobian.setElementAt(20, 18, m31);
5628         jacobian.setElementAt(20, 19, m32);
5629         jacobian.setElementAt(20, 20, m33);
5630 
5631         final var jacobianTrans = jacobian.transposeAndReturnNew();
5632         jacobian.multiply(estimatedCovariance);
5633         jacobian.multiply(jacobianTrans);
5634         estimatedCovariance = jacobian;
5635     }
5636 
5637     /**
5638      * Internal method to perform calibration when common z-axis is assumed for both
5639      * the accelerometer and gyroscope and G-dependent cross biases are ignored.
5640      *
5641      * @throws AlgebraException                              if there are numerical errors.
5642      * @throws FittingException                              if no convergence to solution is found.
5643      * @throws com.irurueta.numerical.NotReadyException      if fitter is not ready.
5644      * @throws InvalidSourceAndDestinationFrameTypeException never happens.
5645      */
5646     private void calibrateCommonAxis() throws AlgebraException, FittingException,
5647             com.irurueta.numerical.NotReadyException, InvalidSourceAndDestinationFrameTypeException {
5648 
5649         // The gyroscope model is
5650         // Ωmeas = bg + (I + Mg) * Ωtrue + Gg * ftrue + w
5651 
5652         // Ideally a least squares solution tries to minimize noise component, so:
5653         // Ωmeas = bg + (I + Mg) * Ωtrue + Gg * ftrue
5654 
5655         // Since G-dependent cross biases are ignored, we can assume that Gg = 0
5656 
5657         // Hence:
5658         // Ωmeas = bg + (I + Mg) * Ωtrue
5659 
5660         // For convergence purposes of the Levenberg-Marquardt algorithm, the
5661         // gyroscope model can be better expressed as:
5662         // Ωmeas = T*K*(Ωtrue + b)
5663         // Ωmeas = M*(Ωtrue + b)
5664         // Ωmeas = M*Ωtrue + M*b
5665 
5666         // where:
5667         // M = I + Mg
5668         // bg = M*b = (I + Mg)*b --> b = M^-1*bg
5669 
5670         // We know that the norm of the true angular rate when the device is in a pixed
5671         // and unknown position and orientation is equal to the Earth rotation rate.
5672         // ||Ωtrue|| = 7.292115E-5 rad/s
5673 
5674         // Hence
5675         // Ωmeas - M*b = M*Ωtrue
5676 
5677         // M^-1 * (Ωmeas - M*b) = Ωtrue
5678 
5679         // ||Ωtrue||^2 = (M^-1 * (Ωmeas - M*b))^T*(M^-1 * (Ωmeas - M*b))
5680         // ||Ωtrue||^2 = (Ωmeas - M*b)^T * (M^-1)^T * M^-1 * (Ωmeas - M*b)
5681         // ||Ωtrue||^2 = (Ωmeas - M*b)^T * ||M^-1||^2 * (Ωmeas - M*b)
5682         // ||Ωtrue||^2 = ||Ωmeas - M*b||^2 * ||M^-1||^2
5683 
5684         // Where:
5685 
5686         // b = [bx]
5687         //     [by]
5688         //     [bz]
5689 
5690         // M = [m11 	m12 	m13]
5691         //     [0 		m22 	m23]
5692         //     [0 	 	0 		m33]
5693 
5694         final var gradientEstimator = new GradientEstimator(this::evaluateCommonAxis);
5695 
5696         final var initialM = Matrix.identity(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
5697         initialM.add(getInitialMg());
5698 
5699         // Force initial M to be upper diagonal
5700         initialM.setElementAt(1, 0, 0.0);
5701         initialM.setElementAt(2, 0, 0.0);
5702         initialM.setElementAt(2, 1, 0.0);
5703 
5704         final var invInitialM = Utils.inverse(initialM);
5705         final var initialBg = getInitialBiasAsMatrix();
5706         final var initialB = invInitialM.multiplyAndReturnNew(initialBg);
5707 
5708         fitter.setFunctionEvaluator(new LevenbergMarquardtMultiDimensionFunctionEvaluator() {
5709             @Override
5710             public int getNumberOfDimensions() {
5711                 // Input points are measured angular rate coordinates
5712                 return BodyKinematics.COMPONENTS;
5713             }
5714 
5715             @Override
5716             public double[] createInitialParametersArray() {
5717                 final var initial = new double[COMMON_Z_AXIS_UNKNOWNS];
5718 
5719                 // biases b
5720                 for (var i = 0; i < BodyKinematics.COMPONENTS; i++) {
5721                     initial[i] = initialB.getElementAtIndex(i);
5722                 }
5723 
5724                 // upper diagonal cross coupling errors M
5725                 var k = BodyKinematics.COMPONENTS;
5726                 for (var j = 0; j < BodyKinematics.COMPONENTS; j++) {
5727                     for (var i = 0; i < BodyKinematics.COMPONENTS; i++) {
5728                         if (i <= j) {
5729                             initial[k] = initialM.getElementAt(i, j);
5730                             k++;
5731                         }
5732                     }
5733                 }
5734 
5735                 return initial;
5736             }
5737 
5738             @Override
5739             public double evaluate(
5740                     final int i, final double[] point, final double[] params, final double[] derivatives)
5741                     throws EvaluationException {
5742 
5743                 measAngularRateX = point[0];
5744                 measAngularRateY = point[1];
5745                 measAngularRateZ = point[2];
5746 
5747                 gradientEstimator.gradient(params, derivatives);
5748 
5749                 return evaluateCommonAxis(params);
5750             }
5751         });
5752 
5753         setInputData();
5754 
5755         fitter.fit();
5756 
5757         final var result = fitter.getA();
5758 
5759         final var bx = result[0];
5760         final var by = result[1];
5761         final var bz = result[2];
5762 
5763         final var m11 = result[3];
5764 
5765         final var m12 = result[4];
5766         final var m22 = result[5];
5767 
5768         final var m13 = result[6];
5769         final var m23 = result[7];
5770         final var m33 = result[8];
5771 
5772         final var mb = new Matrix(BodyKinematics.COMPONENTS, 1);
5773         mb.setElementAtIndex(0, bx);
5774         mb.setElementAtIndex(1, by);
5775         mb.setElementAtIndex(2, bz);
5776 
5777         final var mm = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
5778         mm.setElementAtIndex(0, m11);
5779         mm.setElementAtIndex(1, 0.0);
5780         mm.setElementAtIndex(2, 0.0);
5781 
5782         mm.setElementAtIndex(3, m12);
5783         mm.setElementAtIndex(4, m22);
5784         mm.setElementAtIndex(5, 0.0);
5785 
5786         mm.setElementAtIndex(6, m13);
5787         mm.setElementAtIndex(7, m23);
5788         mm.setElementAtIndex(8, m33);
5789 
5790         setResult(mm, mb);
5791 
5792         // at this point covariance is expressed in terms of b, M and G, and must
5793         // be expressed in terms of bg, Mg and Gg.
5794         // We know that:
5795         // bg = M * b
5796         // Mg = M - I
5797         // Gg = M * G = 0
5798 
5799         // b = [bx]
5800         //     [by]
5801         //     [bz]
5802 
5803         // M = [m11  m12  m13]
5804         //     [0    m22  m23]
5805         //     [0    0    m33]
5806 
5807         // G = [g11  g12  g13] = 0
5808         //     [g21  g22  g23]
5809         //     [g31  g32  g33]
5810 
5811         // bg = [m11  m12  m13][bx] = [m11 * bx + m12 * by + m13 * bz] = [bgx]
5812         //      [0    m22  m23][by]   [           m22 * by + m23 * bz]   [bgy]
5813         //      [0    0    m33][bz]   [                      m33 * bz]   [bgz]
5814 
5815         // Mg = [sx  mxy  mxz] = [m11 - 1    m12         m13     ]
5816         //      [myx sy	  myz]   [0          m22 - 1     m23     ]
5817         //      [mzx mzy  sz ]   [0          0           m33 - 1 ]
5818 
5819         //Gg = [gg11  gg12  gg13] = 0
5820         //     [gg21  gg22  gg23]
5821         //     [gg31  gg32  gg33]
5822 
5823         // Defining the linear application:
5824         // F(b, M) = F(bx, by, bz, m11, m12, m22, m13, m23, m33)
5825         // as:
5826         // [bgx] =  [m11 * bx + m12 * by + m13 * bz]
5827         // [bgy]    [           m22 * by + m23 * bz]
5828         // [bgz]    [                      m33 * bz]
5829         // [sx]     [m11 - 1]
5830         // [sy]     [m22 - 1]
5831         // [sz]     [m33 - 1]
5832         // [mxy]    [m12]
5833         // [mxz]    [m13]
5834         // [myx]    [0]
5835         // [myz]    [m23]
5836         // [mzx]    [0]
5837         // [mzy]    [0]
5838         // [gg11]   [0]
5839         // [gg21]   [0]
5840         // [gg31]   [0]
5841         // [gg12]   [0]
5842         // [gg22]   [0]
5843         // [gg32]   [0]
5844         // [gg13]   [0]
5845         // [gg23]   [0]
5846         // [gg33]   [0]
5847 
5848         // Then the Jacobian of F(b, M) is:
5849         // J = [m11  m12  m13  bx   by   0    bz   0    0  ]
5850         //     [0    m22  m23  0    0    by   0    bz   0  ]
5851         //     [0    0    m33  0    0    0    0    0    bz ]
5852         //     [0    0    0    1    0    0    0    0    0  ]
5853         //     [0    0    0    0    0    1    0    0    0  ]
5854         //     [0    0    0    0    0    0    0    0    1  ]
5855         //     [0    0    0    0    1    0    0    0    0  ]
5856         //     [0    0    0    0    0    0    1    0    0  ]
5857         //     [0    0    0    0    0    0    0    0    0  ]
5858         //     [0    0    0    0    0    0    0    1    0  ]
5859         //     [0    0    0    0    0    0    0    0    0  ]
5860         //     [0    0    0    0    0    0    0    0    0  ]
5861         //     [0    0    0    0    0    0    0    0    0  ]
5862         //     [0    0    0    0    0    0    0    0    0  ]
5863         //     [0    0    0    0    0    0    0    0    0  ]
5864         //     [0    0    0    0    0    0    0    0    0  ]
5865         //     [0    0    0    0    0    0    0    0    0  ]
5866         //     [0    0    0    0    0    0    0    0    0  ]
5867         //     [0    0    0    0    0    0    0    0    0  ]
5868         //     [0    0    0    0    0    0    0    0    0  ]
5869         //     [0    0    0    0    0    0    0    0    0  ]
5870 
5871         // We know that the propagated covariance is J * Cov * J', hence:
5872         final var jacobian = new Matrix(GENERAL_UNKNOWNS_AND_CROSS_BIASES, COMMON_Z_AXIS_UNKNOWNS);
5873 
5874         jacobian.setElementAt(0, 0, m11);
5875         jacobian.setElementAt(0, 1, m12);
5876         jacobian.setElementAt(0, 2, m13);
5877         jacobian.setElementAt(0, 3, bx);
5878         jacobian.setElementAt(0, 4, by);
5879         jacobian.setElementAt(0, 6, bz);
5880 
5881         jacobian.setElementAt(1, 1, m22);
5882         jacobian.setElementAt(1, 2, m23);
5883         jacobian.setElementAt(1, 5, by);
5884         jacobian.setElementAt(1, 7, bz);
5885 
5886         jacobian.setElementAt(2, 2, m33);
5887         jacobian.setElementAt(2, 8, bz);
5888 
5889         jacobian.setElementAt(3, 3, 1.0);
5890         jacobian.setElementAt(4, 5, 1.0);
5891         jacobian.setElementAt(5, 8, 1.0);
5892 
5893         jacobian.setElementAt(6, 4, 1.0);
5894         jacobian.setElementAt(7, 6, 1.0);
5895 
5896         jacobian.setElementAt(9, 7, 1.0);
5897 
5898         final var jacobianTrans = jacobian.transposeAndReturnNew();
5899         jacobian.multiply(estimatedCovariance);
5900         jacobian.multiply(jacobianTrans);
5901         estimatedCovariance = jacobian;
5902     }
5903 
5904     /**
5905      * Internal method to perform general calibration when G-dependent cross biases
5906      * are ignored.
5907      *
5908      * @throws AlgebraException                              if there are numerical errors.
5909      * @throws FittingException                              if no convergence to solution is found.
5910      * @throws com.irurueta.numerical.NotReadyException      if fitter is not ready.
5911      * @throws InvalidSourceAndDestinationFrameTypeException never happens.
5912      */
5913     private void calibrateGeneral() throws AlgebraException, FittingException, com.irurueta.numerical.NotReadyException,
5914             InvalidSourceAndDestinationFrameTypeException {
5915 
5916         // The gyroscope model is
5917         // Ωmeas = bg + (I + Mg) * Ωtrue + Gg * ftrue + w
5918 
5919         // Ideally a least squares solution tries to minimize noise component, so:
5920         // Ωmeas = bg + (I + Mg) * Ωtrue + Gg * ftrue
5921 
5922         // Since G-dependent cross biases are ignored, we can assume that Gg = 0
5923 
5924         // Hence:
5925         // Ωmeas = bg + (I + Mg) * Ωtrue
5926 
5927         // For convergence purposes of the Levenberg-Marquardt algorithm, the
5928         // gyroscope model can be better expressed as:
5929         // Ωmeas = T*K*(Ωtrue + b)
5930         // Ωmeas = M*(Ωtrue + b)
5931         // Ωmeas = M*Ωtrue + M*b
5932 
5933         // where:
5934         // M = I + Mg
5935         // bg = M*b = (I + Mg)*b --> b = M^-1*bg
5936 
5937         // We know that the norm of the true angular rate when the device is in a pixed
5938         // and unknown position and orientation is equal to the Earth rotation rate.
5939         // ||Ωtrue|| = 7.292115E-5 rad/s
5940 
5941         // Hence
5942         // Ωmeas - M*b = M*Ωtrue
5943 
5944         // M^-1 * (Ωmeas - M*b) = Ωtrue
5945 
5946         // ||Ωtrue||^2 = (M^-1 * (Ωmeas - M*b))^T*(M^-1 * (Ωmeas - M*b))
5947         // ||Ωtrue||^2 = (Ωmeas - M*b)^T * (M^-1)^T * M^-1 * (Ωmeas - M*b)
5948         // ||Ωtrue||^2 = (Ωmeas - M*b)^T * ||M^-1||^2 * (Ωmeas - M*b)
5949         // ||Ωtrue||^2 = ||Ωmeas - M*b||^2 * ||M^-1||^2
5950 
5951         // Where:
5952 
5953         // b = [bx]
5954         //     [by]
5955         //     [bz]
5956 
5957         // M = [m11 	m12 	m13]
5958         //     [m21 	m22 	m23]
5959         //     [m31 	m32 	m33]
5960 
5961         final var gradientEstimator = new GradientEstimator(this::evaluateGeneral);
5962 
5963         final var initialM = Matrix.identity(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
5964         initialM.add(getInitialMg());
5965 
5966         final var invInitialM = Utils.inverse(initialM);
5967         final var initialBg = getInitialBiasAsMatrix();
5968         final var initialB = invInitialM.multiplyAndReturnNew(initialBg);
5969 
5970         fitter.setFunctionEvaluator(new LevenbergMarquardtMultiDimensionFunctionEvaluator() {
5971             @Override
5972             public int getNumberOfDimensions() {
5973                 // Input points are measured angular rate coordinates
5974                 return BodyKinematics.COMPONENTS;
5975             }
5976 
5977             @Override
5978             public double[] createInitialParametersArray() {
5979                 final var initial = new double[GENERAL_UNKNOWNS];
5980 
5981                 // biases b
5982                 for (var i = 0; i < BodyKinematics.COMPONENTS; i++) {
5983                     initial[i] = initialB.getElementAtIndex(i);
5984                 }
5985 
5986                 // cross coupling errors M
5987                 final var num = BodyKinematics.COMPONENTS * BodyKinematics.COMPONENTS;
5988                 for (int i = 0, j = BodyKinematics.COMPONENTS; i < num; i++, j++) {
5989                     initial[j] = initialM.getElementAtIndex(i);
5990                 }
5991 
5992                 return initial;
5993             }
5994 
5995             @Override
5996             public double evaluate(
5997                     final int i, final double[] point, final double[] params, final double[] derivatives)
5998                     throws EvaluationException {
5999 
6000                 measAngularRateX = point[0];
6001                 measAngularRateY = point[1];
6002                 measAngularRateZ = point[2];
6003 
6004                 gradientEstimator.gradient(params, derivatives);
6005 
6006                 return evaluateGeneral(params);
6007             }
6008         });
6009 
6010         setInputData();
6011 
6012         fitter.fit();
6013 
6014         final var result = fitter.getA();
6015 
6016         final var bx = result[0];
6017         final var by = result[1];
6018         final var bz = result[2];
6019 
6020         final var m11 = result[3];
6021         final var m21 = result[4];
6022         final var m31 = result[5];
6023 
6024         final var m12 = result[6];
6025         final var m22 = result[7];
6026         final var m32 = result[8];
6027 
6028         final var m13 = result[9];
6029         final var m23 = result[10];
6030         final var m33 = result[11];
6031 
6032         final var mb = new Matrix(BodyKinematics.COMPONENTS, 1);
6033         mb.setElementAtIndex(0, bx);
6034         mb.setElementAtIndex(1, by);
6035         mb.setElementAtIndex(2, bz);
6036 
6037         final var mm = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
6038         mm.setElementAtIndex(0, m11);
6039         mm.setElementAtIndex(1, m21);
6040         mm.setElementAtIndex(2, m31);
6041 
6042         mm.setElementAtIndex(3, m12);
6043         mm.setElementAtIndex(4, m22);
6044         mm.setElementAtIndex(5, m32);
6045 
6046         mm.setElementAtIndex(6, m13);
6047         mm.setElementAtIndex(7, m23);
6048         mm.setElementAtIndex(8, m33);
6049 
6050         setResult(mm, mb);
6051 
6052         // at this point covariance is expressed in terms of b, M and G, and must
6053         // be expressed in terms of bg, Mg and Gg.
6054         // We know that:
6055         // bg = M * b
6056         // Mg = M - I
6057         // Gg = M * G = 0
6058 
6059         // b = [bx]
6060         //     [by]
6061         //     [bz]
6062 
6063         // M = [m11  m12  m13]
6064         //     [m21  m22  m23]
6065         //     [m31  m32  m33]
6066 
6067         // G = [g11  g12  g13] = 0
6068         //     [g21  g22  g23]
6069         //     [g31  g32  g33]
6070 
6071         // bg = [m11  m12  m13][bx] = [m11 * bx + m12 * by + m13 * bz] = [bgx]
6072         //      [m21  m22  m23][by]   [m21 * bx + m22 * by + m23 * bz]   [bgy]
6073         //      [m31  m32  m33][bz]   [m31 * bx + m32 * by + m33 * bz]   [bgz]
6074 
6075         // Mg = [sx  mxy  mxz] = [m11 - 1    m12         m13     ]
6076         //      [myx sy	  myz]   [m21        m22 - 1     m23     ]
6077         //      [mzx mzy  sz ]   [m31        m32         m33 - 1 ]
6078 
6079         // Gg = [gg11  gg12  gg13] = 0
6080         //      [gg21  gg22  gg23]
6081         //      [gg31  gg32  gg33]
6082 
6083         // Defining the linear application:
6084         // F(b, M) = F(bx, by, bz, m11, m21, m31, m12, m22, m32, m13, m23, m33)
6085         // as:
6086         // [bgx] =  [m11 * bx + m12 * by + m13 * bz]
6087         // [bgy]    [m21 * bx + m22 * by + m23 * bz]
6088         // [bgz]    [m31 * bx + m32 * by + m33 * bz]
6089         // [sx]     [m11 - 1]
6090         // [sy]     [m22 - 1]
6091         // [sz]     [m33 - 1]
6092         // [mxy]    [m12]
6093         // [mxz]    [m13]
6094         // [myx]    [m21]
6095         // [myz]    [m23]
6096         // [mzx]    [m31]
6097         // [mzy]    [m32]
6098         // [gg11]   [0]
6099         // [gg21]   [0]
6100         // [gg31]   [0]
6101         // [gg12]   [0]
6102         // [gg22]   [0]
6103         // [gg32]   [0]
6104         // [gg13]   [0]
6105         // [gg23]   [0]
6106         // [gg33]   [0]
6107 
6108         // Then the Jacobian of F(b, M) is:
6109         // J = [m11  m12  m13  bx   0    0    by   0    0    bz   0    0 ]
6110         //     [m21  m22  m23  0    bx   0    0    by   0    0    bz   0 ]
6111         //     [m31  m32  m33  0    0    bx   0    0    by   0    0    bz]
6112         //     [0    0    0    1    0    0    0    0    0    0    0    0 ]
6113         //     [0    0    0    0    0    0    0    1    0    0    0    0 ]
6114         //     [0    0    0    0    0    0    0    0    0    0    0    1 ]
6115         //     [0    0    0    0    0    0    1    0    0    0    0    0 ]
6116         //     [0    0    0    0    0    0    0    0    0    1    0    0 ]
6117         //     [0    0    0    0    1    0    0    0    0    0    0    0 ]
6118         //     [0    0    0    0    0    0    0    0    0    0    1    0 ]
6119         //     [0    0    0    0    0    1    0    0    0    0    0    0 ]
6120         //     [0    0    0    0    0    0    0    0    1    0    0    0 ]
6121         //     [0    0    0    0    0    0    0    0    0    0    0    0 ]
6122         //     [0    0    0    0    0    0    0    0    0    0    0    0 ]
6123         //     [0    0    0    0    0    0    0    0    0    0    0    0 ]
6124         //     [0    0    0    0    0    0    0    0    0    0    0    0 ]
6125         //     [0    0    0    0    0    0    0    0    0    0    0    0 ]
6126         //     [0    0    0    0    0    0    0    0    0    0    0    0 ]
6127         //     [0    0    0    0    0    0    0    0    0    0    0    0 ]
6128         //     [0    0    0    0    0    0    0    0    0    0    0    0 ]
6129         //     [0    0    0    0    0    0    0    0    0    0    0    0 ]
6130 
6131         // We know that the propagated covariance is J * Cov * J', hence:
6132         final var jacobian = new Matrix(GENERAL_UNKNOWNS_AND_CROSS_BIASES, GENERAL_UNKNOWNS);
6133 
6134         jacobian.setElementAt(0, 0, m11);
6135         jacobian.setElementAt(0, 1, m12);
6136         jacobian.setElementAt(0, 2, m13);
6137         jacobian.setElementAt(0, 3, bx);
6138         jacobian.setElementAt(0, 6, by);
6139         jacobian.setElementAt(0, 9, bz);
6140 
6141         jacobian.setElementAt(1, 0, m21);
6142         jacobian.setElementAt(1, 1, m22);
6143         jacobian.setElementAt(1, 2, m23);
6144         jacobian.setElementAt(1, 4, bx);
6145         jacobian.setElementAt(1, 7, by);
6146         jacobian.setElementAt(1, 10, bz);
6147 
6148         jacobian.setElementAt(2, 0, m31);
6149         jacobian.setElementAt(2, 1, m32);
6150         jacobian.setElementAt(2, 2, m33);
6151         jacobian.setElementAt(2, 5, bx);
6152         jacobian.setElementAt(2, 8, by);
6153         jacobian.setElementAt(2, 11, bz);
6154 
6155         jacobian.setElementAt(3, 3, 1.0);
6156         jacobian.setElementAt(4, 7, 1.0);
6157         jacobian.setElementAt(5, 11, 1.0);
6158 
6159         jacobian.setElementAt(6, 6, 1.0);
6160         jacobian.setElementAt(7, 9, 1.0);
6161         jacobian.setElementAt(8, 4, 1.0);
6162 
6163         jacobian.setElementAt(9, 10, 1.0);
6164         jacobian.setElementAt(10, 5, 1.0);
6165         jacobian.setElementAt(11, 8, 1.9);
6166 
6167         final var jacobianTrans = jacobian.transposeAndReturnNew();
6168         jacobian.multiply(estimatedCovariance);
6169         jacobian.multiply(jacobianTrans);
6170         estimatedCovariance = jacobian;
6171     }
6172 
6173     /**
6174      * Sets input data into Levenberg-Marquardt fitter when G-dependent cross biases
6175      * are taken into account.
6176      *
6177      * @throws AlgebraException                              if provided accelerometer cross coupling
6178      *                                                       errors are not valid.
6179      * @throws InvalidSourceAndDestinationFrameTypeException never happens
6180      */
6181     private void setInputDataWithGDependentCrossBiases() throws AlgebraException,
6182             InvalidSourceAndDestinationFrameTypeException {
6183         // compute reference frame at current position
6184         final var nedPosition = getNedPosition();
6185         final var nedC = new CoordinateTransformation(FrameType.BODY_FRAME, FrameType.LOCAL_NAVIGATION_FRAME);
6186         final var nedFrame = new NEDFrame(nedPosition, nedC);
6187         final var ecefFrame = NEDtoECEFFrameConverter.convertNEDtoECEFAndReturnNew(nedFrame);
6188         final var refKinematics = ECEFKinematicsEstimator.estimateKinematicsAndReturnNew(timeInterval, ecefFrame,
6189                 ecefFrame);
6190 
6191         final var refAngularRateX = refKinematics.getAngularRateX();
6192         final var refAngularRateY = refKinematics.getAngularRateY();
6193         final var refAngularRateZ = refKinematics.getAngularRateZ();
6194 
6195         final var w2 = turntableRotationRate * turntableRotationRate;
6196 
6197         final var numMeasurements = measurements.size();
6198         final var x = new Matrix(numMeasurements, 2 * BodyKinematics.COMPONENTS);
6199         final var y = new double[numMeasurements];
6200         final var angularRateStandardDeviations = new double[numMeasurements];
6201         var i = 0;
6202         for (final var measurement : measurements) {
6203             final var measuredKinematics = measurement.getKinematics();
6204 
6205             final var angularRateX = measuredKinematics.getAngularRateX();
6206             final var angularRateY = measuredKinematics.getAngularRateY();
6207             final var angularRateZ = measuredKinematics.getAngularRateZ();
6208 
6209             final var fX = measuredKinematics.getFx();
6210             final var fY = measuredKinematics.getFy();
6211             final var fZ = measuredKinematics.getFz();
6212 
6213             x.setElementAt(i, 0, angularRateX - refAngularRateX);
6214             x.setElementAt(i, 1, angularRateY - refAngularRateY);
6215             x.setElementAt(i, 2, angularRateZ - refAngularRateZ);
6216 
6217             x.setElementAt(i, 3, fX);
6218             x.setElementAt(i, 4, fY);
6219             x.setElementAt(i, 5, fZ);
6220 
6221             y[i] = w2;
6222 
6223             angularRateStandardDeviations[i] = measurement.getAngularRateStandardDeviation();
6224 
6225             i++;
6226         }
6227 
6228         fitter.setInputData(x, y, angularRateStandardDeviations);
6229 
6230         ba = getAccelerometerBiasAsMatrix();
6231         ma = getAccelerometerMa();
6232         accelerationFixer.setBias(ba);
6233         accelerationFixer.setCrossCouplingErrors(ma);
6234     }
6235 
6236     /**
6237      * Sets input data into Levenberg-Marquardt fitter when G-dependent cross biases
6238      * are ignored.
6239      *
6240      * @throws AlgebraException                              if provided accelerometer cross coupling
6241      *                                                       errors are not valid.
6242      * @throws InvalidSourceAndDestinationFrameTypeException never happens.
6243      */
6244     private void setInputData() throws AlgebraException, InvalidSourceAndDestinationFrameTypeException {
6245 
6246         // compute reference frame at current position
6247         final var nedPosition = getNedPosition();
6248         final var nedC = new CoordinateTransformation(FrameType.BODY_FRAME, FrameType.LOCAL_NAVIGATION_FRAME);
6249         final var nedFrame = new NEDFrame(nedPosition, nedC);
6250         final var ecefFrame = NEDtoECEFFrameConverter.convertNEDtoECEFAndReturnNew(nedFrame);
6251         final var refKinematics = ECEFKinematicsEstimator.estimateKinematicsAndReturnNew(timeInterval, ecefFrame,
6252                 ecefFrame);
6253 
6254         final var refAngularRateX = refKinematics.getAngularRateX();
6255         final var refAngularRateY = refKinematics.getAngularRateY();
6256         final var refAngularRateZ = refKinematics.getAngularRateZ();
6257 
6258         final var w2 = turntableRotationRate * turntableRotationRate;
6259 
6260         final var numMeasurements = measurements.size();
6261         final var x = new Matrix(numMeasurements, BodyKinematics.COMPONENTS);
6262         final var y = new double[numMeasurements];
6263         final var angularRateStandardDeviations = new double[numMeasurements];
6264         var i = 0;
6265         for (final var measurement : measurements) {
6266             final var measuredKinematics = measurement.getKinematics();
6267 
6268             final var angularRateX = measuredKinematics.getAngularRateX();
6269             final var angularRateY = measuredKinematics.getAngularRateY();
6270             final var angularRateZ = measuredKinematics.getAngularRateZ();
6271 
6272             x.setElementAt(i, 0, angularRateX - refAngularRateX);
6273             x.setElementAt(i, 1, angularRateY - refAngularRateY);
6274             x.setElementAt(i, 2, angularRateZ - refAngularRateZ);
6275 
6276             y[i] = w2;
6277 
6278             angularRateStandardDeviations[i] = measurement.getAngularRateStandardDeviation();
6279 
6280             i++;
6281         }
6282 
6283         fitter.setInputData(x, y, angularRateStandardDeviations);
6284 
6285         ba = getAccelerometerBiasAsMatrix();
6286         ma = getAccelerometerMa();
6287         accelerationFixer.setBias(ba);
6288         accelerationFixer.setCrossCouplingErrors(ma);
6289     }
6290 
6291     /**
6292      * Converts provided NED position expressed in terms of latitude, longitude and height respect
6293      * mean Earth surface, to position expressed in ECEF coordinates.
6294      *
6295      * @param position NED position to be converted.
6296      * @return converted position expressed in ECEF coordinates.
6297      */
6298     private static ECEFPosition convertPosition(final NEDPosition position) {
6299         final var velocity = new ECEFVelocity();
6300         final var result = new ECEFPosition();
6301         NEDtoECEFPositionVelocityConverter.convertNEDtoECEF(
6302                 position.getLatitude(), position.getLongitude(), position.getHeight(), 0.0, 0.0, 0.0,
6303                 result, velocity);
6304         return result;
6305     }
6306 
6307     /**
6308      * Converts acceleration instance to meters per squared second.
6309      *
6310      * @param acceleration acceleration instance to be converted.
6311      * @return converted value.
6312      */
6313     private static double convertAcceleration(final Acceleration acceleration) {
6314         return AccelerationConverter.convert(acceleration.getValue().doubleValue(), acceleration.getUnit(),
6315                 AccelerationUnit.METERS_PER_SQUARED_SECOND);
6316     }
6317 
6318     /**
6319      * Converts angular speed value and unit to radians per second.
6320      *
6321      * @param value angular speed value.
6322      * @param unit  unit of angular speed value.
6323      * @return converted value.
6324      */
6325     private static double convertAngularSpeed(final double value, final AngularSpeedUnit unit) {
6326         return AngularSpeedConverter.convert(value, unit, AngularSpeedUnit.RADIANS_PER_SECOND);
6327     }
6328 
6329     /**
6330      * Converts angular speed instance to radians per second.
6331      *
6332      * @param angularSpeed angular speed instance to be converted.
6333      * @return converted value.
6334      */
6335     private static double convertAngularSpeed(final AngularSpeed angularSpeed) {
6336         return convertAngularSpeed(angularSpeed.getValue().doubleValue(), angularSpeed.getUnit());
6337     }
6338 
6339     /**
6340      * Converts time instance to seconds.
6341      *
6342      * @param time time instance to be converted.
6343      * @return converted value.
6344      */
6345     private static double convertTime(final Time time) {
6346         return TimeConverter.convert(time.getValue().doubleValue(), time.getUnit(), TimeUnit.SECOND);
6347     }
6348 
6349     /**
6350      * Makes proper conversion of internal cross-coupling, bias and g-dependent
6351      * cross bias matrices.
6352      *
6353      * @param m internal scaling and cross-coupling matrix.
6354      * @param b internal bias matrix.
6355      * @param g internal g-dependent cross bias matrix.
6356      * @throws AlgebraException if a numerical instability occurs.
6357      */
6358     private void setResult(final Matrix m, final Matrix b, final Matrix g) throws AlgebraException {
6359         setResult(m, b);
6360 
6361         // Gg = M*G
6362         m.multiply(g, estimatedGg);
6363     }
6364 
6365     /**
6366      * Makes proper conversion of internal cross-coupling and bias matrices.
6367      *
6368      * @param m internal scaling and cross-coupling matrix.
6369      * @param b internal bias matrix.
6370      * @throws AlgebraException if a numerical instability occurs.
6371      */
6372     private void setResult(final Matrix m, final Matrix b) throws AlgebraException {
6373         // Because:
6374         // M = I + Mg
6375         // b = M^-1*bg
6376 
6377         // Then:
6378         // Mg = M - I
6379         // bg = M*b
6380 
6381         if (estimatedBiases == null) {
6382             estimatedBiases = new double[BodyKinematics.COMPONENTS];
6383         }
6384 
6385         final var bg = m.multiplyAndReturnNew(b);
6386         bg.toArray(estimatedBiases);
6387 
6388         if (estimatedMg == null) {
6389             estimatedMg = m;
6390         } else {
6391             estimatedMg.copyFrom(m);
6392         }
6393 
6394         for (var i = 0; i < BodyKinematics.COMPONENTS; i++) {
6395             estimatedMg.setElementAt(i, i, estimatedMg.getElementAt(i, i) - 1.0);
6396         }
6397 
6398         if (estimatedGg == null) {
6399             estimatedGg = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
6400         } else {
6401             estimatedGg.initialize(0.0);
6402         }
6403 
6404         estimatedCovariance = fitter.getCovar();
6405         estimatedChiSq = fitter.getChisq();
6406         estimatedChiSqDegreesOfFreedom = fitter.getChisqDegreesOfFreedom();
6407         estimatedReducedChiSq = estimatedChiSq / (double) estimatedChiSqDegreesOfFreedom;
6408         estimatedMse = fitter.getMse();
6409         try {
6410             estimatedP = fitter.getP();
6411             estimatedQ = fitter.getQ();
6412         } catch (final MaxIterationsExceededException | IllegalArgumentException ignore) {
6413             // only happens for numerically unstable values or if the number of estimated parameters is larger than the
6414             // number of samples
6415             estimatedP = 1.0;
6416             estimatedQ = 0.0;
6417         }
6418     }
6419 
6420     /**
6421      * Computes estimated true angular rate squared norm using current measured
6422      * angular rate and specific force along with provided parameters for the
6423      * general case when G-dependent cross biases are taken into account.
6424      * This method is internally executed during gradient estimation and
6425      * Levenberg-Marquardt fitting needed for calibration computation.
6426      *
6427      * @param params array containing parameters for the general purpose case
6428      *               when G-dependent cross biases are taken into account. Must
6429      *               have length 21.
6430      * @return estimated true angular rate squared norm.
6431      * @throws EvaluationException if there are numerical instabilities.
6432      */
6433     private double evaluateGeneralWithGDependentCrossBiases(final double[] params) throws EvaluationException {
6434         final var bx = params[0];
6435         final var by = params[1];
6436         final var bz = params[2];
6437 
6438         final var m11 = params[3];
6439         final var m21 = params[4];
6440         final var m31 = params[5];
6441 
6442         final var m12 = params[6];
6443         final var m22 = params[7];
6444         final var m32 = params[8];
6445 
6446         final var m13 = params[9];
6447         final var m23 = params[10];
6448         final var m33 = params[11];
6449 
6450         final var g11 = params[12];
6451         final var g21 = params[13];
6452         final var g31 = params[14];
6453 
6454         final var g12 = params[15];
6455         final var g22 = params[16];
6456         final var g32 = params[17];
6457 
6458         final var g13 = params[18];
6459         final var g23 = params[19];
6460         final var g33 = params[20];
6461 
6462         return evaluate(bx, by, bz, m11, m21, m31, m12, m22, m32, m13, m23, m33,
6463                 g11, g21, g31, g12, g22, g32, g13, g23, g33);
6464     }
6465 
6466     /**
6467      * Computes estimated true angular rate squared norm using current measured
6468      * angular rate and specific force along with provided parameters when
6469      * common z-axis is assumed and G-dependent cross biases are taken into
6470      * account.
6471      * This method is internally executed during gradient estimation and
6472      * Levenberg-Marquardt fitting needed for calibration computation.
6473      *
6474      * @param params array containing parameters for the general purpose case
6475      *               when G-dependent cross biases are taken into account. Must
6476      *               have length 18.
6477      * @return estimated true angular rate squared norm.
6478      * @throws EvaluationException if there are numerical instabilities.
6479      */
6480     private double evaluateCommonAxisWithGDependentCrossBiases(final double[] params) throws EvaluationException {
6481         final var bx = params[0];
6482         final var by = params[1];
6483         final var bz = params[2];
6484 
6485         final var m11 = params[3];
6486 
6487         final var m12 = params[4];
6488         final var m22 = params[5];
6489 
6490         final var m13 = params[6];
6491         final var m23 = params[7];
6492         final var m33 = params[8];
6493 
6494         final var g11 = params[9];
6495         final var g21 = params[10];
6496         final var g31 = params[11];
6497 
6498         final var g12 = params[12];
6499         final var g22 = params[13];
6500         final var g32 = params[14];
6501 
6502         final var g13 = params[15];
6503         final var g23 = params[16];
6504         final var g33 = params[17];
6505 
6506         return evaluate(bx, by, bz, m11, 0.0, 0.0, m12, m22, 0.0, m13, m23, m33,
6507                 g11, g21, g31, g12, g22, g32, g13, g23, g33);
6508     }
6509 
6510     /**
6511      * Computes estimated true angular rate squared norm using current measured
6512      * angular rate and provided parameters for the general case when G-dependent
6513      * cross biases are ignored.
6514      * This method is internally executed during gradient estimation and
6515      * Levenberg-Marquardt fitting needed for calibration computation.
6516      *
6517      * @param params array containing current parameters for the general purpose case
6518      *               when G-dependent cross biases are ignored. Must have length 12.
6519      * @return estimated true angular rate squared norm.
6520      * @throws EvaluationException if there are numerical instabilities.
6521      */
6522     private double evaluateGeneral(final double[] params) throws EvaluationException {
6523         final var bx = params[0];
6524         final var by = params[1];
6525         final var bz = params[2];
6526 
6527         final var m11 = params[3];
6528         final var m21 = params[4];
6529         final var m31 = params[5];
6530 
6531         final var m12 = params[6];
6532         final var m22 = params[7];
6533         final var m32 = params[8];
6534 
6535         final var m13 = params[9];
6536         final var m23 = params[10];
6537         final var m33 = params[11];
6538 
6539         return evaluate(bx, by, bz, m11, m21, m31, m12, m22, m32, m13, m23, m33);
6540     }
6541 
6542     /**
6543      * Computes estimated true angular rate squared norm using current measured
6544      * angular rate and provided parameters when common z-axis is assumed and
6545      * G-dependent cross biases are ignored.
6546      * This method is internally executed during gradient estimation and
6547      * Levenberg-Marquardt fitting needed for calibration computation.
6548      *
6549      * @param params array containing current parameters for the common z-axis case
6550      *               when G-dependent cross biases are ignored. Must have length 9.
6551      * @return estimated true angular rate squared norm.
6552      * @throws EvaluationException if there are numerical instabilities.
6553      */
6554     private double evaluateCommonAxis(final double[] params) throws EvaluationException {
6555         final var bx = params[0];
6556         final var by = params[1];
6557         final var bz = params[2];
6558 
6559         final var m11 = params[3];
6560 
6561         final var m12 = params[4];
6562         final var m22 = params[5];
6563 
6564         final var m13 = params[6];
6565         final var m23 = params[7];
6566         final var m33 = params[8];
6567 
6568         return evaluate(bx, by, bz, m11, 0.0, 0.0, m12, m22, 0.0, m13, m23, m33);
6569     }
6570 
6571     /**
6572      * Computes estimated true angular rate squared norm using current measured
6573      * angular rate and provided parameters.
6574      * This method is internally executed during gradient estimation and
6575      * Levenberg-Marquardt fitting needed for calibration computation.
6576      *
6577      * @param bx  x-coordinate of bias.
6578      * @param by  y-coordinate of bias.
6579      * @param bz  z-coordinate of bias.
6580      * @param m11 element 1,1 of cross-coupling error matrix.
6581      * @param m21 element 2,1 of cross-coupling error matrix.
6582      * @param m31 element 3,1 of cross-coupling error matrix.
6583      * @param m12 element 1,2 of cross-coupling error matrix.
6584      * @param m22 element 2,2 of cross-coupling error matrix.
6585      * @param m32 element 3,2 of cross-coupling error matrix.
6586      * @param m13 element 1,3 of cross-coupling error matrix.
6587      * @param m23 element 2,3 of cross-coupling error matrix.
6588      * @param m33 element 3,3 of cross-coupling error matrix.
6589      * @param g11 element 1,1 of g-dependent cross bias matrix.
6590      * @param g21 element 2,1 of g-dependent cross bias matrix.
6591      * @param g31 element 3,1 of g-dependent cross bias matrix.
6592      * @param g12 element 1,2 of g-dependent cross bias matrix.
6593      * @param g22 element 2,2 of g-dependent cross bias matrix.
6594      * @param g32 element 3,2 of g-dependent cross bias matrix.
6595      * @param g13 element 1,3 of g-dependent cross bias matrix.
6596      * @param g23 element 2,3 of g-dependent cross bias matrix.
6597      * @param g33 element 3,3 of g-dependent cross bias matrix.
6598      * @return estimated true angular rate squared norm.
6599      * @throws EvaluationException if there are numerical instabilities.
6600      */
6601     private double evaluate(final double bx, final double by, final double bz,
6602                             final double m11, final double m21, final double m31,
6603                             final double m12, final double m22, final double m32,
6604                             final double m13, final double m23, final double m33,
6605                             final double g11, final double g21, final double g31,
6606                             final double g12, final double g22, final double g32,
6607                             final double g13, final double g23, final double g33) throws EvaluationException {
6608 
6609         // Ωmeas = bg + (I + Mg) * Ωtrue + Gg * ftrue
6610         // Ωmeas = M*(Ωtrue + b + G * ftrue)
6611 
6612         // M = I + Mg
6613         // bg = M*b --> b = M^-1*bg
6614         // Gg = M*G --> G = M^-1*Gg
6615 
6616         // Ωtrue = M^-1 * Ωmeas - b - G*ftrue
6617 
6618         try {
6619             if (measAngularRate == null) {
6620                 measAngularRate = new Matrix(BodyKinematics.COMPONENTS, 1);
6621             }
6622             if (fmeas == null) {
6623                 fmeas = new Matrix(BodyKinematics.COMPONENTS, 1);
6624             }
6625             if (m == null) {
6626                 m = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
6627             }
6628             if (invM == null) {
6629                 invM = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
6630             }
6631             if (b == null) {
6632                 b = new Matrix(BodyKinematics.COMPONENTS, 1);
6633             }
6634             if (g == null) {
6635                 g = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
6636             }
6637             if (trueAngularRate == null) {
6638                 trueAngularRate = new Matrix(BodyKinematics.COMPONENTS, 1);
6639             }
6640             if (ftrue == null) {
6641                 ftrue = new Matrix(BodyKinematics.COMPONENTS, 1);
6642             }
6643             if (ba == null) {
6644                 ba = new Matrix(BodyKinematics.COMPONENTS, 1);
6645             }
6646             if (ma == null) {
6647                 ma = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
6648             }
6649             if (tmp == null) {
6650                 tmp = new Matrix(BodyKinematics.COMPONENTS, 1);
6651             }
6652 
6653             measAngularRate.setElementAtIndex(0, measAngularRateX);
6654             measAngularRate.setElementAtIndex(1, measAngularRateY);
6655             measAngularRate.setElementAtIndex(2, measAngularRateZ);
6656 
6657             fmeas.setElementAtIndex(0, fmeasX);
6658             fmeas.setElementAtIndex(1, fmeasY);
6659             fmeas.setElementAtIndex(2, fmeasZ);
6660 
6661             m.setElementAt(0, 0, m11);
6662             m.setElementAt(1, 0, m21);
6663             m.setElementAt(2, 0, m31);
6664 
6665             m.setElementAt(0, 1, m12);
6666             m.setElementAt(1, 1, m22);
6667             m.setElementAt(2, 1, m32);
6668 
6669             m.setElementAt(0, 2, m13);
6670             m.setElementAt(1, 2, m23);
6671             m.setElementAt(2, 2, m33);
6672 
6673             Utils.inverse(m, invM);
6674 
6675             b.setElementAtIndex(0, bx);
6676             b.setElementAtIndex(1, by);
6677             b.setElementAtIndex(2, bz);
6678 
6679             g.setElementAt(0, 0, g11);
6680             g.setElementAt(1, 0, g21);
6681             g.setElementAt(2, 0, g31);
6682 
6683             g.setElementAt(0, 1, g12);
6684             g.setElementAt(1, 1, g22);
6685             g.setElementAt(2, 1, g32);
6686 
6687             g.setElementAt(0, 2, g13);
6688             g.setElementAt(1, 2, g23);
6689             g.setElementAt(2, 2, g33);
6690 
6691             getAccelerometerBiasAsMatrix(ba);
6692             getAccelerometerMa(ma);
6693 
6694             // fix measured accelerometer value to obtain true
6695             // specific force
6696             accelerationFixer.fix(fmeas, ftrue);
6697             g.multiply(ftrue, tmp);
6698 
6699             invM.multiply(measAngularRate, trueAngularRate);
6700             trueAngularRate.subtract(b);
6701             trueAngularRate.subtract(tmp);
6702 
6703             final var norm = Utils.normF(trueAngularRate);
6704             return norm * norm;
6705 
6706         } catch (final AlgebraException e) {
6707             throw new EvaluationException(e);
6708         }
6709     }
6710 
6711     /**
6712      * Computes estimated true angular rate squared norm using current measured
6713      * angular rate and provided parameters.
6714      * This method is internally executed during gradient estimation and
6715      * Levenberg-Marquardt fitting needed for calibration computation.
6716      *
6717      * @param bx  x-coordinate of bias.
6718      * @param by  y-coordinate of bias.
6719      * @param bz  z-coordinate of bias.
6720      * @param m11 element 1,1 of cross-coupling error matrix.
6721      * @param m21 element 2,1 of cross-coupling error matrix.
6722      * @param m31 element 3,1 of cross-coupling error matrix.
6723      * @param m12 element 1,2 of cross-coupling error matrix.
6724      * @param m22 element 2,2 of cross-coupling error matrix.
6725      * @param m32 element 3,2 of cross-coupling error matrix.
6726      * @param m13 element 1,3 of cross-coupling error matrix.
6727      * @param m23 element 2,3 of cross-coupling error matrix.
6728      * @param m33 element 3,3 of cross-coupling error matrix.
6729      * @return estimated true angular rate squared norm.
6730      * @throws EvaluationException if there are numerical instabilities.
6731      */
6732     private double evaluate(final double bx, final double by, final double bz,
6733                             final double m11, final double m21, final double m31,
6734                             final double m12, final double m22, final double m32,
6735                             final double m13, final double m23, final double m33) throws EvaluationException {
6736 
6737         // Ωmeas = M*(Ωtrue + b)
6738         // Ωtrue = M^-1 * Ωmeas - b
6739 
6740         try {
6741             if (measAngularRate == null) {
6742                 measAngularRate = new Matrix(BodyKinematics.COMPONENTS, 1);
6743             }
6744             if (m == null) {
6745                 m = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
6746             }
6747             if (invM == null) {
6748                 invM = new Matrix(BodyKinematics.COMPONENTS, BodyKinematics.COMPONENTS);
6749             }
6750             if (b == null) {
6751                 b = new Matrix(BodyKinematics.COMPONENTS, 1);
6752             }
6753             if (trueAngularRate == null) {
6754                 trueAngularRate = new Matrix(BodyKinematics.COMPONENTS, 1);
6755             }
6756 
6757             measAngularRate.setElementAtIndex(0, measAngularRateX);
6758             measAngularRate.setElementAtIndex(1, measAngularRateY);
6759             measAngularRate.setElementAtIndex(2, measAngularRateZ);
6760 
6761             m.setElementAt(0, 0, m11);
6762             m.setElementAt(1, 0, m21);
6763             m.setElementAt(2, 0, m31);
6764 
6765             m.setElementAt(0, 1, m12);
6766             m.setElementAt(1, 1, m22);
6767             m.setElementAt(2, 1, m32);
6768 
6769             m.setElementAt(0, 2, m13);
6770             m.setElementAt(1, 2, m23);
6771             m.setElementAt(2, 2, m33);
6772 
6773             Utils.inverse(m, invM);
6774 
6775             b.setElementAtIndex(0, bx);
6776             b.setElementAtIndex(1, by);
6777             b.setElementAtIndex(2, bz);
6778 
6779             invM.multiply(measAngularRate, trueAngularRate);
6780             trueAngularRate.subtract(b);
6781 
6782             final var norm = Utils.normF(trueAngularRate);
6783             return norm * norm;
6784 
6785         } catch (final AlgebraException e) {
6786             throw new EvaluationException(e);
6787         }
6788     }
6789 }