Kalman estimation of position and velocity for ReaTHM testing applications

Offshore wind power research is a rapidly growing field, because of the present climate crisis and increasing focus on renewable energy. Model testing plays an important role in the risk and cost analysis associated with offshore wind turbines (OWTs). The real-time hybrid model testing concept (ReaTHM testing) solves important challenges related to model testing of OWTs, such as achieving an accurate modelling of the wind field, and the occurrence of scaling issues when modelling wind and waves simultaneously. However, ReaTHM test set-ups are generally sensitive to noise, signal loss and inaccuracies in sensor values. The present study is focused on the design and implementation of a state estimator able to accurately estimate the position and velocity of floating structures, while taking disturbances into account. By combining the information received from several different sensors with mathematical models, the estimator provides smooth and reliable position and velocity estimates for ReaTHM testing applications. The main objective of the present study is to develop a kinematic state space model that could represent the motion of any floating structure in six degrees of freedom (6-DOF). The kinematic model is implemented in MATLAB, and acceleration time series obtained with numerical simulations are used as inputs. The computed outputs agree with the corresponding simulated motions. A Kalman estimator based on the state space model is designed, implemented and tested against virtual data from the numerical model, with artificially added disturbances. Sensitivity analyses addressing the robustness towards noise, time delays, signal loss and uncertainties are performed to identify the limits of the estimator. The estimator is demonstrated to be robust to most types of disturbances. Further, the state estimator is tested against physical data from laboratory experiments. Good agreement between the physically measured and the estimated states is observed.


Introduction
The recent energy and climate crisis has led to an increasing focus towards the renewable energy, and offshore wind is one of the most important renewable energy sources. Compared to onshore wind power projects, offshore wind turbines (OWTs) require significant extra engineering and research efforts on installation, component design, electrical connection and materials [1]. On the other hand, they provide the opportunity for developments in areas previously not IOP Conf. Series: Journal of Physics: Conf. Series 1104 (2018) 012008 IOP Publishing doi: 10.1088/1742-6596/1104/1/012008 2 Figure 1. Real-time hybrid model testing approach deemed suitable for wind turbines, giving access to a high energy potential without disturbing neighbouring communities. The production and operational costs of OWTs has been reduced in recent years due to the introduction of bigger turbines and higher voltage cables [2]. However, the costs need to be further reduced, compared with other renewable energy technologies, for OWTs to be commercially viable. Limited knowledge on the environmental loads and the structural responses of wind turbines induces risks and uncertainties, which are the major cost driving factors in the design of OWTs. Model testing can play an important role in decreasing the impact of these factors, and is also a way of detecting over-conservative designs (thus reducing costs). However, there are several challenges related to the model testing of OWTs [3]. These include inherent modelling difficulties, issues linked to the simultaneous physical modelling of wind and waves, and difficulties in implementing the control systems for blade angles and measurements.
A possible solution to these challenges is to use real-time hybrid model (ReaTHM) testing, as described by Sauder et al [4]. In ReaTHM testing, the system under study is divided into a physical and a numerical substructure that interact via sensors and actuators integrated in the physical substructure (Fig. 1). This means that wind and aerodynamic loads may be replaced by actuators that are controlled by outputs from a numerical model. The inputs to the numerical model are the real-time measured motions of the structure, which are used to compute the rotor loads in full scale. These loads are then scaled down to model scale, and applied to the structure in real time. In this way, the Froude-scaled mass properties of the whole turbine are maintained. The ReaTHM testing concept is a subclass of hardware-in-the-loop testing [4]. Such test set-ups are generally sensitive to noise, signal loss and inaccuracies in sensor values. Any erroneous input values to the numerical substructure may elicit erroneous or harmful actuator control signals. The focus of the present work is therefore to design a state estimator that estimates and filters the positions and velocities of the physical substructure, to achieve better accuracy and reliability than what is possible through measurements alone.
A method for ReaTHM testing of OWTs is proposed by Chabaud [5], who suggests a design and verification procedure for ReaTHM test campaigns and provides simple position and velocity observer designs. Although these observers worked satisfactorily in that particular study, an optimal state estimator (Kalman estimator) would typically be better suited for suppressing erroneous properties in positions and velocities. It would use the past errors between predictions and measurements to provide optimal estimates, hence minimizing the estimation error in the statistical sense. In the marine control system designed by Vilsen et al. [6], a non-linear passive observer (NPO) was designed and tested. The objectives of those tests were to generate estimates of the states, and to filter out noise from the signals.
In the current work, values and variations in attitude angles (i.e. roll, pitch, yaw) are assumed to be small, so that a linear estimation approach can be used. This enables the use of Kalman estimators, which rely on a linear system model. The estimation performance is tested by artificially adding different types of disturbances. The estimators are further validated using physical data from previous experiments, conducted by Vilsen et al. [6].

Numerical Model
In this section, a high-level overview of the system design is given as well as more detailed descriptions of the kinematic modelling and the estimator design.

Overview
Two different versions of the system are designed for tests using virtual and physical data, respectively (Fig. 2). Virtual data provided by the simulation platform SIMA [8] is used for the development and verification of the estimator, while physical data obtained from laboratory experiments is required for the final validation. The simulator that is meant to represent the physical system does not have exactly the same inputs and outputs as the physical system. Therefore, two sets of system matrices are necessary.
A kinematic state-space model is designed to represent the motion of any floating structure in 6-DOF. The virtual acceleration time series (U ) is given as input and the corresponding system states and outputs (X and Y , respectively) are computed. It is desirable to maintain these time series without any disturbances, for comparison purposes. A plant model that takes process noise originating from resonant structural vibrations (v 1 ) into account is therefore implemented using the same state-space matrices. The plant model is intended to simulate the physical system, and measurement noise (v 2 ) is therefore added to the output of the plant model, yielding the measured outputs plantY . A Kalman estimator with U and plantY as inputs is then designed using the same state-space matrices. The intent of the Kalman estimator is to provide estimates of the states and outputs (kalmX and kalmY , respectively). Fig. 2a shows the overall system design for tests using virtual data, where acceleration (U ) and output from the plant model (plantY ) are used as inputs to the Kalman estimator.
For the tests using physical data, it is not necessary to artificially add noise or other disturbances, since these are inherently present in the physical data. The plant model is therefore not needed, meaning that the kinematic model is only used to derive the system matrices of the Kalman estimator. The physical acceleration time series (U ) is given as input to the Kalman estimator, along with the process outputs trueY , which are subject to real measurement and process noise. The overall system design for tests using physical data is shown in Fig. 2b.

Kinematic model
Since the estimator should be applicable to any floating structure, the mathematical models incorporated into the estimator should be limited to the kinematic relations between positions, velocities and accelerations. The state vector consists of the variables to be estimated: 4 where x, y, and z are the translation coordinates of the structure's center of mass (COM), and φ, θ, and ψ are the Euler angles that define the structure's orientation. Both translation and rotation are defined with respect to the inertial frame (NED), and the order of rotations is yaw(ψ), pitch (θ) and then roll (φ). The Euler angles are chosen to represent the attitude angles, because the angles are assumed small and singularities will not occur. The variables u, v, and w correspond to the linear velocities of the COM, while p, q and r denote the angular velocities of the COM. For a time instant t, the Euler angles and the angular velocities of the model are related through the equation: where s · = sin(·), c · = cos(·) and t · = tan(·). Angular velocities and linear accelerations are measured through gyrometers and accelerometers attached to the structure, and are hence defined with respect to BODY. The translations and orientations are provided by an optical position measurement system (OPMS), and are therefore defined with respect to NED. The variables which can be measured, i.e. the translations, orientations, angular velocities and linear accelerations (a x , a y , and a z ) of the physical substructure, are considered as outputs. The output vector is therefore: Therefore, the following continuous-time state-space system is proposed: where v 1 and v 2 are process and output noise, respectively. Since the angular accelerationsṗ,q, andṙ are not measured, they are considered as process noise. The angular accelerations should only add to the differentiated angular velocities, hence v 1 = [0 9×1 ,ṗ,q,ṙ] T . The input vector is defined as U = [a x , a y , a z ] T , while the matrices A(X), B, C, D, G, and H are defined according to Fossen [7], section 2.2.1, yielding: The non-linearity of the system through the entries R n b (Θ nb ) and T Θ (Θ nb ) in the matrix A(X) should be noticed. However, the attitude angles involved (i.e. roll, pitch, yaw) will normally be small, so the linear approximations described by Fossen [7], section 2.2.1, could be used.

Simplified model for testing with SIMA
The state estimator is first tested against a SIMA model of a floating wind turbine developed by Chabaud [5]. The SIMA software simulates the aerodynamic and hydrodynamic loads (wind, waves and current), and computes the global position, velocity and acceleration data of the structure [8]. 5 Some simplifications and modifications of the kinematic model are necessary for this case study. Both the linear and the angular accelerations are now treated as inputs to the kinematic model, since they can be accessed in the SIMA output data. This yields the input vector U = [a x , a y , a z ,ṗ,q,ṙ] T . The initial condition of the state vector is taken as the first set of position and velocity values from SIMA. All positions, velocities and accelerations are defined globally in SIMA, unlike in the physical tests, where linear accelerations and angular velocities are measured locally. In the simplified version, all values are therefore defined in the same (global) reference frame. Hence, the rotation and transformation matrices R n b (θ nb ) and T θ (θ nb ) are replaced by identity matrices (I). This makes the simplified kinematic model linear and time-invariant. The continuous-time state-space system is then defined as follows: The kinematic model is implemented in MATLAB, and acceleration time series obtained with the SIMA simulations are used as inputs. The computed outputs agree with the corresponding simulated motions (results not shown here, but can be found in [9], section 6.1), hence providing a verification of the kinematic model defined in 2.2.

Estimator design
The kinematic model forms the basis of the estimator design. There are many applicationspecific approaches for estimating an unknown state from a set of process measurements, but many of these methods do not consider the typically noisy nature of the measurements [10]. In the present study, the fundamental sources of information are always the same: position, angular velocity and linear acceleration measurements. These are all derived from noisy sensors, and the noise is typically statistical in nature (or can be effectively modelled as such). This implies that stochastic methods should be used to address the problem. Since the angular motion is relatively small, the system is approximated as linear. The solution adopted in this study is therefore the Kalman estimator, which would be more appropriate for ReaTHM testing applications than the state observers previously implemented [5,6] since it uses the previous errors between predictions and measurements to provide optimal estimates. The state-space model is defined in continuous time because this allows the use of the kinematic relations between position, velocity and acceleration. However, the Kalman estimator must in practice be implemented in discrete time. It is therefore decided to discretize the kinematic model and perform all computations in discrete time. The next states and outputs are found using the discrete system dynamics: In the Kalman estimator design, the process and measurement noise vectors are assumed to be white, Gaussian and zero-mean: The time-varying one-step ahead Kalman predictor is given by [11]: where is a gain matrix "weighting" the innovation, e(k) = C d v(k) + v 2 (k). The state prediction error can be expressed as: and the state prediction error covariance matrix can be computed recursively with the Ricatti equation: A convenient simplification when there is no expected change in the noise covariance matrix is to replace K(k) and P (k + 1) with their steady-state values: The constant matricesK andP can then be computed off-line beforehand. This yields the steady-state one-step ahead Kalman predictor: This is a sub-optimal estimator, but it allows the practical handling of high-dimensional problems. The steady-state gain matrix can be computed as: while the solution of the Algebraic Riccati Equation (ARE) yields the matrixP : Both steady-state and time-varying versions of the Kalman estimator are designed, implemented in MATLAB and tested. Time-varying Kalman estimators perform better than steady-state Kalman estimators in some cases, since they take changes in the noise covariance into account.

Results and discussion
This section presents the functionality of the Kalman estimator quantitatively. The numerical model is validated, first using virtual data and then physical data. Due to limited space, only the key results will be shown here. The full test data sets are available in the Github repository [12] and are discussed in details by Mehammer [9].  The sensitivity of the estimator to different types of noise is investigated, and the results for measurement noise is shown here to illustrate (Fig. 3). Measurement noise can only be seen at the outputs, not at the states. Fig. 3a presents the time evolution of the first output from the kinematic model (actual), the first output from the plant model (measured), and the first output from the Kalman estimator (estimated) when 10% white measurement noise is added. The estimator is able to track the actual output in the presence of measurement noise. A high peak is observed in both signals around 20 s, then the values damp out with time. The relative estimation error (absolute estimation error divided by absolute value) for state 1 (x-translation) in the presence of 10% measurement noise (Fig. 3b) stays below 0.15% throughout most of the simulation, but has two high peaks at around t=0 s and t=280 s. This is because the xtranslation has low absolute values at these time points. The estimator is thereby shown to be robust towards measurement noise.  Further, tests are carried out to investigate the sensitivity of the estimator to different kinds of uncertainties. As an illustration, Fig. 4 presents the sensitivity to uncertain values of the measurement noise covariance matrix R, in the presence of 1% measurement noise and 1% process noise. In Fig. 4a, the time evolutions of the actual, measured and estimated output 1 (x-translation) is shown for an R matrix with values ten times larger values than the ideal noise information. It can be observed that uncertain measurement noise values have limited effect on the estimation performance, since the estimated output 1 is smooth and accurate. The relative  Fig. 4b. A peak can be observed at the beginning of the simulation. After this, the relative difference between the actual and the estimated state is insignificant until t=270 s. Then it increases with time, reaching a peak value of around 4.5% at t=470 s. Although recent improvements have reduced the time delays between the optical and inertial measurements in ReaHTM tests to beneath 20 ms, delays were typically in the range 20-40 ms at the time of the present study. Hence, simulations to test the effects of time delays are carried out with a delay of 1 time step (50 ms). Fig. 5a presents the actual, measured and estimated values of state 1 with a delay of 1 time step. Good agreement is observed between all three signals. A few high peaks can be observed in the relative estimation error (Fig. 5b) due to the low absolute values of the x-translation at these time values. Otherwise, the relative estimation error stays below 0.5%. Hence, the estimator is shown to be robust towards time delays. Lastly, the sensitivity to signal loss is investigated. The measured position outputs are set to zero to simulate loss of the OPMS signal. Fig. 6a presents the time evolutions of the first output of the kinematic model (actual), the first output of the plant model (measured) and the first output of the Kalman estimator (estimated) when the position measurements are lost for 10 consecutive time steps at four different occasions (t=100 s, t=200 s, t=300 s, t=400 s). The measured output drops to zero at these values of time, while the estimated output stays close to the actual output. Hence, the estimator is robust to signal loss of 0.5 s duration for this output (x-translation). Fig. 6b presents the relative estimation error for output 1 in the presence of signal loss of duration 10 time steps. It has peaks at around t=0 s and t=280 s, due to the low absolute values of the x-translation at these values of time.

Validation of the Kalman estimator using physical data
Further, the state estimator is tested against physical data from the laboratory experiments conducted by Vilsen et al. [6]. The experimental data was obtained using a buoy that was 9 moored to the walls of a still-water basin at three points. The buoy was pulled away from its equilibrium position and released, leading to a decaying oscillatory movement as the buoy gradually converged back to the equilibrium due to the stiffness of the mooring system. Both the steady-state and the time-varying implementations of the Kalman estimator are tested, and knowledge about delays and inaccuracies in the sensors used is taken into account. The OPMS delay is compensated for by predicting the position data 25 ms ahead of time. Furthermore, the impact of gravity is compensated for by rotating the measured accelerations to the global reference frame, and then subtracting the gravitational acceleration in the z-direction. The validation is performed by comparing the estimator results with the experimental data. Good results are obtained for both versions of the Kalman estimator, as implied by the results for state 1 and state 7 (linear velocity in x-direction).   Fig. 7a presents the time evolution of the measured state 1, along with the estimations computed by the steady-state and the time-varying Kalman estimators (only a short section of the signals are shown, in order to highlight the interesting parts). Good agreement is observed between all three signals. As expected, the relative estimation error is relatively high when the measured value of the x-translation is close to zero (Fig. 7b). However, the error stays below 2% for the parts of the measured signal with relevant magnitude. The relative estimation errors are similar for the steady-state and the time-varying Kalman estimators. The linear velocities are not measured in the physical tests. To evaluate the performance of the two different versions of the Kalman estimator, the linear velocities are therefore compared with the values estimated by the non-linear passive observer (NPO) implemented by Vilsen et al. [6]. Fig. 8a  10 values decay uniformly. When magnifying a short section of the graph, it is apparent that the Kalman estimates are smoother than the NPO estimates, and that the NPO has a longer time delay than the Kalman estimators (Fig. 8b). This might be because the Kalman estimators provide optimal estimates, minimizing the estimation error in the statistical sense.

Conclusions and further work
The current project focuses on developing a kinematic state estimator using a sensor fusion algorithm, and a Kalman estimator is found to provide reliable position and velocity estimates in 6-DOF for ReaTHM testing applications. The main task of the estimator is to handle disturbances such as noise, signal loss, time delays and uncertainties. Sensitivity analyses are carried out to identify the limits of the Kalman estimator. Virtual data from a floating wind turbine case study is used in the sensitivity analyses. The estimator is further validated through comparison with physical data. The main conclusions of the work are: • The generic kinematic model developed can recreate the SIMA simulated motions with reasonable accuracy. • A Kalman estimator providing smooth and accurate position and velocity estimates in six degrees of freedom is designed, implemented and tested. • The Kalman estimator is proven to be robust towards disturbances such as measurement noise, uncertainties, time delays and signal loss. • The Kalman estimator is able to estimate the states with a good accuracy, when compared with physical measurements. An improvement from the previously implemented estimators is demonstrated.
In the further work, ReaTHM tests should be carried out to investigate the estimator performance in closed-loop and real-time applications. In such tests, Kalman filter tuning should be performed to adjust the noise covariance matrices Q, R and N , according to the different applications. For cases where the angular motion is more significant, an extended Kalman filter (EKF) or an NPO should be implemented to better deal with the non-linearities.