# Discrete Wavelet Transform-based Approach of Real-time Wave Filtering for Dynamic Positioning of Marine Vessels

Authors: Changde Liu; Yufang Zhang; Min Gu; Longhui Zhang; Yanbin Teng; Zhengfeng Liu; Qiang Sun; Naxin Wei
Journal: Mathematical Problems in Engineering
 (2022)

---

## Abstract

The ongoing development of deep-sea resources has contributed toward the widespread use of dynamic positioning (DP) systems that can operate in arbitrary sea areas without the limitation of operating water depth compared to mooring systems. The second-order slowly varying forces can be compensated by DP feedback control. Wave filtering techniques are crucial for designing the DP controller, which can improve the positioning performance and service of multiple actuators. The objective of wave filtering is to separate the motion induced by the first-order high-frequency oscillatory waves from the motion caused by the second-order low-frequency (LF) slowly varying disturbances. This paper proposes a novel real-time wave filtering approach for DP systems on the basis of the discrete wavelet transform without relying on the previous knowledge of dynamics of the vessel. Meanwhile, to suppress the end effects of wavelet filtering and achieve real-time application, symmetric boundary extension and the sliding window method are investigated in order to capture the slow time-varying LF motion of DP vessels. The proposed real-time wavelet filtering approach can significantly improve the quality of the separated components by selecting the optimal decomposition level, wavelet function parameters, and threshold denoising method. The real-time and filtering performance of the proposed approach is verified using a scaled model of an offshore supply vessel as well as full-scale experiments on the EXPLOERE I scientific investigation vessel under the action of winds, currents, and waves. The experimental results confirm that the proposed wave filtering algorithm is efficient and that the LF component waveform can be effectively recovered after wavelet filtering.


---

## Body

## 1. Introduction

Dynamic positioning (DP) systems are defined as a set of components for maintaining a floating vessel in a specific position or along a predefined track by means of sufficient thrusters and propellers to ensure that the position and heading of the vessel are unaffected by the disturbances to which it is subjected. Compared to the traditional vessel positioning operation mode using mooring systems, DP systems can perform well in deep and shallow water, and its costs do not increase significantly with the water depth. Moreover, DP systems can arrive at the designated sea area more rapidly and effectively and perform the required operation immediately. Meanwhile, the higher redundancy for ensuring the safety of DP systems has resulted in higher positioning accuracy and increased payload compared to mooring systems. At present, DP systems are used not only in the oil industry but also in various marine applications such as drilling, pipe laying, and supply.For surface marine vessels with DP systems, only the horizontal plan position and heading are controlled. The environmental disturbances mainly include winds, waves, and ocean currents. In a DP system, only the low-frequency (LF) estimation of the vessel motion must be counteracted by the thrust system; the rapid and purely oscillatory wave-frequency (WF) motion (first-order wave-induced load) should not enter the feedback control system because the WF components in the applied thrust may endanger the propulsion system by causing excessive wear, resulting in shorter propulsion unit service intervals and service life. When designing DP systems for vessels, a key issue to consider is wave filtering. In most cases, the position and heading measurements are corrupted with sensor noise.In the first DP systems invented for marine vessels in the 1960s, the conventional PID controller with a low-pass filter was employed to remove the first-order wave-induced motion [1]. However, low-pass filters introduce phase lag, and the filtered data are processed backward using a low-pass filter to reduce the phase lag [2]. Subsequently, the performance of wave filtering was improved using a Kalman filter [3–6]. The DP model proposed in [4] has the advantages of self-adaptive noise filtering, optimum combination of different measures of the position and heading of the vessel, and effective estimates of the environmental forces. However, this approach is based on a linear mathematical model of the dynamic system, which is only valid for small heading changes (around 10°); hence, the lineation process of the state-space model must be carried out at different heading angles between 0° and 360°. To ensure effective DP control for different operation conditions, the matrix parameters must be tuned with the change in the heading angle. Thus, to enhance the robustness of wave filtering for different heading angles, the extended Kalman filter (EKF) was developed [7, 8], and it remains in use in existing DP commercial applications, such as the Kongsberg DP system. Compared with the KF and EKF, H∞ filtering has the advantage of robustness against the unmodeled dynamics. A state-space H∞ filter was proposed in [9]. Recently, many studies were reviewed in [10]. The H∞ filter has been applied to DP, and it has a similar potential compared with EKF [11]. In the late 20th century, with the increasing application of observer techniques, a nonlinear observer with wave filtering and estimation was designed using passivity methods [12]. In [13], this observer was extended to adaptive wave filtering. In particular, second-order wave transfer function approximations have been used extensively [14–18]. A new linear design model was employed for DP vessels with wave filtering by suggesting that the WF motions should be modeled in a hydrodynamic frame. Meanwhile, a second-order wave transfer function approximation and a modified model for the WF components of motion were proposed; the gains and design of the observer were based on the WF model [19, 20], and three different environmental conditions were considered. There is a need for an observer that achieves wave filtering and separates the LF and WF position and heading. However, all these observer designs assume that the DP system operates in a stable sea state, or a priori wave model parameters have been established by designing separate observers from calm to high seas [21]. The aforementioned assumption is not realistic, as the sea state is constantly changing, and the mass gain-scheduled parameters must be tuned to satisfy different sea states. Therefore, for the ideal observer design, it is difficult to automatically adjust the reconstruction in accordance with the varying sea state.To mitigate the aforementioned difficulties, an adaptive filtering approach using empirical mode decomposition (EMD) was adopted for wave filtering of a DP vessel, and experimental results were presented in [22]. The obvious advantage of the EMD filter approach is that it is model-independent, as it requires no information of the system dynamics. Although the EMD method can remove the WF motion components, useful LF motion components may be lost because of the overlapping of different mode functions. The signals from a full-scale oil tank were analyzed using an offline filtering algorithm, and online EMD of wave filtering was proposed [22]. In [23], the performance of the online EMD filter with a sliding window was evaluated via DP numerical simulation, and it was proposed that the time window and sliding window lengths should be reduced for real-time DP applications.At present, few studies on wavelet filtering have been conducted in the field of vessel motion signal processing [24, 25]. Wavelet filtering has two main advantages. First, it is possible to filter short signals. Second, there is no significant signal delay caused by the filtration process. In a previous study, we verified that the wavelet filtering method can filter the measured noise to improve the precision of ship motion prediction under irregular waves [24]. Wavelet threshold denoising of signals has also been adopted in the identification of the maneuvering model on the basis of preprocessing data corrupted by high-frequency measured noise in calm water [25].In view of the aforementioned advantages of the wavelet filtering algorithm, this paper first proposes a new wave filtering approach for DP vessels based on the discrete wavelet transform (DWT), which is model-independent and only needs information of the wave characteristics in advance compared to the traditional commercial DP system using EKF. In the actual operation process of the DP system under different sea states, the acquisition of the wave characteristics is much easier compared to dynamic modeling, and the modeling difficulties caused by different vessel types under hypothetical conditions are avoided. Next, to suppress the end effects of wavelet filtering and achieve real-time application, symmetric boundary extension and the sliding window method are investigated in order to capture the slow time-varying LF motion of DP vessels. Another advantage of the proposed method is that it is hardly affected by phase lag and amplitude distortion owing to the use of near-symmetric symlet wavelets. Finally, the effectiveness and feasibility of the proposed wave filtering approach are verified using a scaled model of an offshore supply vessel as well as experiments on a full-scale scientific investigation vessel under the action of winds, currents, and waves. The experimental results confirm that the wave filtering algorithm is suitable for analyzing nonstationary and nonlinear signals when the position and heading of the DP vessel contain high-frequency components.The remainder of this paper is organized as follows. Section2 introduces the mathematical problem considered in this study, including the wavelet transform, decomposition, and threshold filtering steps. Section 3 describes the scaled and full-scale experimental set-up as well as the experimental environment conditions. Section 4 presents, discusses, and analyzes the experimental results of the scaled and full-scale vessels. Finally, Section 5 concludes the paper.


## 2. Methods

### 2.1. Wavelet Transform

A wavelet is essentially a small wave that grows and decays over a limited period. In the continuous wavelet transform (CWT), a given signal of finite energy is projected to a continuous family of frequency bands or subspaces of various scales inL2R. The wavelet transform is visualized in the continuous domain as follows:(1)WTfa,b=f,ψa,b=∫−∞∞ftψ¯a,btdt,where ft∈L2R is an arbitrary function to be transformed, and ψ¯a,bt is the conjugate function of ψa,bt.For instance, the signal may be represented in each subspace of a scale for alla. For each a, the subspace is generated by the translation and rescaling of a single function ψt∈L2R called the mother wavelet which is given as follows:(2)ψa,bt=a−1/2ψt−ba,where a,b∈R and a≠0; a denotes the scaling parameter and b denotes the translation parameter. Further, a represents the time and frequency resolutions of the scaled mother wavelet, while b represents its translation along the time axis.The continuous scale is computationally intractable; thus, DWT comes into play. DWT minimizes the transformation for discrete values ofa and b while guaranteeing invertibility [26]. To use this transform for sampling data, the scaling and shifting parameters are discretized on a sampling grid as follows:(3)a=a0m,b=nb0a0m,m,n∈Z,a0>1,b0>0.Thus, we get(4)ψm,nt=a0−m/2ψt−nb0a0ma0m=a0−m/2ψa0−mt−nb0.Notably, the shifting parameterb is a function of the scale a. This represents a crucial advantage of DWT, which corresponds to the sampling of a and b such that the consecutive discrete values of a and b as well as the sampling intervals differ by a factor of two. More formally, a0=2 and b0=1. The discrete wavelet is used as shown in the equation below, which constitutes a family of orthonormal basis functions:(5)ψm,nt=2−m/2ψ2−mt−n.In general, the wavelet coefficients can be derived accordingly as follows:(6)cm,nf=ft,ψm,nt=2−m/2∫−∞∞ftψ2−mt−ndt.Further,ft can be reconstructed using the discrete wavelet coefficients as follows:(7)ft=∑m∑ncm,nfψm,ntdt.When the input functionft and the wavelet parameters a and b are represented in a discrete form, the transformation is referred to as the DWT of the signal [27].


### 2.2. Wavelet Decomposition

DWT extracts the time and frequency localized coefficients from an input signal or function, which constitute its wavelet decomposition. Mallet proposed a fast algorithm to compute an input signal by using the wavelet multiresolution algorithm [27]. This algorithm is based on orthogonal multiresolution analysis.According to the scalability and inclusiveness of multiresolution analysis, wavelet decompositions are conducted from wavelet basis functions that consist of a scaling functionϕt and a wavelet function ψt. These wavelet decompositions decompose a signal into a series of low- and high-frequency sub-bands. An orthonormal basis is formed by the scaling and wavelet functions in the Hilbert space (H) [28, 29]:(8)ϕm,nt=2−m/2ϕ2−mt−n,ψm,nt=2−m/2ψ2−mt−n.m,n∈Z,ψt∈L2R,Any functionft∈H can be represented as(9)ft=∑k=−∞∞cj0,kϕj0,kt+∑j=j0∞∑k=−∞∞dj,kψj,kt,where j0 denotes the initial decomposition level, and cj,k and dj,k denote the approximation and detail coefficients at scale j and location k, respectively.


### 2.3. Threshold Processing of Coefficients

A typical threshold processing function consists of hard and soft value processing functions [30].For hard thresholding, it can be defined as follows:(10)y=x,x≥σ,0,x<σ.For soft thresholding, it can be defined as follows:(11)y=sgnxgx−σ,x≥σ,0,x<σ.In general,g is a linear function, and σ denotes the threshold value.


### 2.4. Real-Time Wavelet Filter

#### 2.4.1. Wavelet Filtering

During wavelet decomposition, the signal is decomposed into a series of low- and high-frequency sub-bands via scaling and wavelet functions. Successive lower-frequency band-pass filters (wavelets) are used to orthogonally decompose the signal, and the lowest band is analyzed using the scaling function. The coefficients in this domain are removed such that only the local modulus maxima are retained. The signal is then reconstructed with the new coefficients to produce a filtering version of the original signal.Realizing the signal filter using the wavelet transform involves three basic steps as follows:Step 1: Wavelet transform for decomposing the data in the trend part, applied to the original time series data to obtain the same number of coefficients as the size of the data.Step 2: Determination of the optimal decomposition level.Step 3: Threshold filtering method of approximation and detail coefficients. The new coefficients are used to regenerate the signal using the inverse DWT (IDWT).Decomposition and reconstruction can also be performed using Mallet’s pyramid algorithm [27]. Figure 1 shows an example of the overall process of decomposition and reconstruction by the wavelet transform.Figure 1
Process of decomposition and reconstruction by the wavelet transform.In Figure1, cj,k represents the coarse approximation (LF trend), and dj,k represents the detailed information (HF noise). The difference between the first-level approximation coefficients c1,k and the original series x yields the detail coefficients of the first-level d1,k. To obtain c2,k, c1,k is approximated by a set of basis functions. The coefficients that pass through the low-pass filter are called approximation coefficients, and the coefficients that pass through the high-pass filter are called detail coefficients.


#### 2.4.2. Symmetric Boundary Extension

DWT is implemented using a two-band perfect reconstruction filter bank. In many applications, such as signal filtering [31, 32], compactly supported wavelets are required to realize perfect signal reconstruction. However, a real orthogonal symmetric wavelet basis, such as the Haar wavelet, is not continuous and the filter is of the order 1 only, which is not adequate for practical applications. The design of the ideal orthogonal and symmetric wavelet filter is challenging. For applications of signal processing in DP experiments, we can select orthogonal and near-symmetric wavelets to have a minimal phase response rather than an exact linear phase. Therefore, an orthogonal and near-symmetric wavelet of the symlet wavelet family is selected in this study. The effect of the optimal decomposition level and wavelet family will be discussed in Section 4.In practical engineering applications, iterative decomposition and reconstruction are required when using the mallet algorithm. Therefore, the end effects of the filtered signal are inevitable because the input sequence of the signal has a finite length [33]. To suppress the end effects of wavelet filtering for real-time application, symmetric boundary extension (SBE) is applied to wavelet filtering in this study. The extension procedure is as follows. Assuming that xnn=0,1,…,N−1 is the N point vector of the input samples, for the sequence xn, the input signal is extended by adding N/2 (N-1/2) values when N is even (odd) at each end, and then, the new sequences x^n can be constructed by SBE as follows:(12)x^n=x−n−1,for−N/2≤n<0,x2N−n−1,forN≤n<3N/2.


#### 2.4.3. Sliding Window Method

To realize the real-time application of wavelet filtering, a sliding window method [34] is adopted to capture the slow time-varying LF motion of DP vessels. The development of the sliding window is shown in Figure 2. The sliding window is fixed, and the real-time filtered result is obtained using the newly measured and historical data.Figure 2
Schematic diagram of the sliding window.


## 2.1. Wavelet Transform

A wavelet is essentially a small wave that grows and decays over a limited period. In the continuous wavelet transform (CWT), a given signal of finite energy is projected to a continuous family of frequency bands or subspaces of various scales inL2R. The wavelet transform is visualized in the continuous domain as follows:(1)WTfa,b=f,ψa,b=∫−∞∞ftψ¯a,btdt,where ft∈L2R is an arbitrary function to be transformed, and ψ¯a,bt is the conjugate function of ψa,bt.For instance, the signal may be represented in each subspace of a scale for alla. For each a, the subspace is generated by the translation and rescaling of a single function ψt∈L2R called the mother wavelet which is given as follows:(2)ψa,bt=a−1/2ψt−ba,where a,b∈R and a≠0; a denotes the scaling parameter and b denotes the translation parameter. Further, a represents the time and frequency resolutions of the scaled mother wavelet, while b represents its translation along the time axis.The continuous scale is computationally intractable; thus, DWT comes into play. DWT minimizes the transformation for discrete values ofa and b while guaranteeing invertibility [26]. To use this transform for sampling data, the scaling and shifting parameters are discretized on a sampling grid as follows:(3)a=a0m,b=nb0a0m,m,n∈Z,a0>1,b0>0.Thus, we get(4)ψm,nt=a0−m/2ψt−nb0a0ma0m=a0−m/2ψa0−mt−nb0.Notably, the shifting parameterb is a function of the scale a. This represents a crucial advantage of DWT, which corresponds to the sampling of a and b such that the consecutive discrete values of a and b as well as the sampling intervals differ by a factor of two. More formally, a0=2 and b0=1. The discrete wavelet is used as shown in the equation below, which constitutes a family of orthonormal basis functions:(5)ψm,nt=2−m/2ψ2−mt−n.In general, the wavelet coefficients can be derived accordingly as follows:(6)cm,nf=ft,ψm,nt=2−m/2∫−∞∞ftψ2−mt−ndt.Further,ft can be reconstructed using the discrete wavelet coefficients as follows:(7)ft=∑m∑ncm,nfψm,ntdt.When the input functionft and the wavelet parameters a and b are represented in a discrete form, the transformation is referred to as the DWT of the signal [27].


## 2.2. Wavelet Decomposition

DWT extracts the time and frequency localized coefficients from an input signal or function, which constitute its wavelet decomposition. Mallet proposed a fast algorithm to compute an input signal by using the wavelet multiresolution algorithm [27]. This algorithm is based on orthogonal multiresolution analysis.According to the scalability and inclusiveness of multiresolution analysis, wavelet decompositions are conducted from wavelet basis functions that consist of a scaling functionϕt and a wavelet function ψt. These wavelet decompositions decompose a signal into a series of low- and high-frequency sub-bands. An orthonormal basis is formed by the scaling and wavelet functions in the Hilbert space (H) [28, 29]:(8)ϕm,nt=2−m/2ϕ2−mt−n,ψm,nt=2−m/2ψ2−mt−n.m,n∈Z,ψt∈L2R,Any functionft∈H can be represented as(9)ft=∑k=−∞∞cj0,kϕj0,kt+∑j=j0∞∑k=−∞∞dj,kψj,kt,where j0 denotes the initial decomposition level, and cj,k and dj,k denote the approximation and detail coefficients at scale j and location k, respectively.


## 2.3. Threshold Processing of Coefficients

A typical threshold processing function consists of hard and soft value processing functions [30].For hard thresholding, it can be defined as follows:(10)y=x,x≥σ,0,x<σ.For soft thresholding, it can be defined as follows:(11)y=sgnxgx−σ,x≥σ,0,x<σ.In general,g is a linear function, and σ denotes the threshold value.


## 2.4. Real-Time Wavelet Filter

### 2.4.1. Wavelet Filtering

During wavelet decomposition, the signal is decomposed into a series of low- and high-frequency sub-bands via scaling and wavelet functions. Successive lower-frequency band-pass filters (wavelets) are used to orthogonally decompose the signal, and the lowest band is analyzed using the scaling function. The coefficients in this domain are removed such that only the local modulus maxima are retained. The signal is then reconstructed with the new coefficients to produce a filtering version of the original signal.Realizing the signal filter using the wavelet transform involves three basic steps as follows:Step 1: Wavelet transform for decomposing the data in the trend part, applied to the original time series data to obtain the same number of coefficients as the size of the data.Step 2: Determination of the optimal decomposition level.Step 3: Threshold filtering method of approximation and detail coefficients. The new coefficients are used to regenerate the signal using the inverse DWT (IDWT).Decomposition and reconstruction can also be performed using Mallet’s pyramid algorithm [27]. Figure 1 shows an example of the overall process of decomposition and reconstruction by the wavelet transform.Figure 1
Process of decomposition and reconstruction by the wavelet transform.In Figure1, cj,k represents the coarse approximation (LF trend), and dj,k represents the detailed information (HF noise). The difference between the first-level approximation coefficients c1,k and the original series x yields the detail coefficients of the first-level d1,k. To obtain c2,k, c1,k is approximated by a set of basis functions. The coefficients that pass through the low-pass filter are called approximation coefficients, and the coefficients that pass through the high-pass filter are called detail coefficients.


### 2.4.2. Symmetric Boundary Extension

DWT is implemented using a two-band perfect reconstruction filter bank. In many applications, such as signal filtering [31, 32], compactly supported wavelets are required to realize perfect signal reconstruction. However, a real orthogonal symmetric wavelet basis, such as the Haar wavelet, is not continuous and the filter is of the order 1 only, which is not adequate for practical applications. The design of the ideal orthogonal and symmetric wavelet filter is challenging. For applications of signal processing in DP experiments, we can select orthogonal and near-symmetric wavelets to have a minimal phase response rather than an exact linear phase. Therefore, an orthogonal and near-symmetric wavelet of the symlet wavelet family is selected in this study. The effect of the optimal decomposition level and wavelet family will be discussed in Section 4.In practical engineering applications, iterative decomposition and reconstruction are required when using the mallet algorithm. Therefore, the end effects of the filtered signal are inevitable because the input sequence of the signal has a finite length [33]. To suppress the end effects of wavelet filtering for real-time application, symmetric boundary extension (SBE) is applied to wavelet filtering in this study. The extension procedure is as follows. Assuming that xnn=0,1,…,N−1 is the N point vector of the input samples, for the sequence xn, the input signal is extended by adding N/2 (N-1/2) values when N is even (odd) at each end, and then, the new sequences x^n can be constructed by SBE as follows:(12)x^n=x−n−1,for−N/2≤n<0,x2N−n−1,forN≤n<3N/2.


### 2.4.3. Sliding Window Method

To realize the real-time application of wavelet filtering, a sliding window method [34] is adopted to capture the slow time-varying LF motion of DP vessels. The development of the sliding window is shown in Figure 2. The sliding window is fixed, and the real-time filtered result is obtained using the newly measured and historical data.Figure 2
Schematic diagram of the sliding window.


## 2.4.1. Wavelet Filtering

During wavelet decomposition, the signal is decomposed into a series of low- and high-frequency sub-bands via scaling and wavelet functions. Successive lower-frequency band-pass filters (wavelets) are used to orthogonally decompose the signal, and the lowest band is analyzed using the scaling function. The coefficients in this domain are removed such that only the local modulus maxima are retained. The signal is then reconstructed with the new coefficients to produce a filtering version of the original signal.Realizing the signal filter using the wavelet transform involves three basic steps as follows:Step 1: Wavelet transform for decomposing the data in the trend part, applied to the original time series data to obtain the same number of coefficients as the size of the data.Step 2: Determination of the optimal decomposition level.Step 3: Threshold filtering method of approximation and detail coefficients. The new coefficients are used to regenerate the signal using the inverse DWT (IDWT).Decomposition and reconstruction can also be performed using Mallet’s pyramid algorithm [27]. Figure 1 shows an example of the overall process of decomposition and reconstruction by the wavelet transform.Figure 1
Process of decomposition and reconstruction by the wavelet transform.In Figure1, cj,k represents the coarse approximation (LF trend), and dj,k represents the detailed information (HF noise). The difference between the first-level approximation coefficients c1,k and the original series x yields the detail coefficients of the first-level d1,k. To obtain c2,k, c1,k is approximated by a set of basis functions. The coefficients that pass through the low-pass filter are called approximation coefficients, and the coefficients that pass through the high-pass filter are called detail coefficients.


## 2.4.2. Symmetric Boundary Extension

DWT is implemented using a two-band perfect reconstruction filter bank. In many applications, such as signal filtering [31, 32], compactly supported wavelets are required to realize perfect signal reconstruction. However, a real orthogonal symmetric wavelet basis, such as the Haar wavelet, is not continuous and the filter is of the order 1 only, which is not adequate for practical applications. The design of the ideal orthogonal and symmetric wavelet filter is challenging. For applications of signal processing in DP experiments, we can select orthogonal and near-symmetric wavelets to have a minimal phase response rather than an exact linear phase. Therefore, an orthogonal and near-symmetric wavelet of the symlet wavelet family is selected in this study. The effect of the optimal decomposition level and wavelet family will be discussed in Section 4.In practical engineering applications, iterative decomposition and reconstruction are required when using the mallet algorithm. Therefore, the end effects of the filtered signal are inevitable because the input sequence of the signal has a finite length [33]. To suppress the end effects of wavelet filtering for real-time application, symmetric boundary extension (SBE) is applied to wavelet filtering in this study. The extension procedure is as follows. Assuming that xnn=0,1,…,N−1 is the N point vector of the input samples, for the sequence xn, the input signal is extended by adding N/2 (N-1/2) values when N is even (odd) at each end, and then, the new sequences x^n can be constructed by SBE as follows:(12)x^n=x−n−1,for−N/2≤n<0,x2N−n−1,forN≤n<3N/2.


## 2.4.3. Sliding Window Method

To realize the real-time application of wavelet filtering, a sliding window method [34] is adopted to capture the slow time-varying LF motion of DP vessels. The development of the sliding window is shown in Figure 2. The sliding window is fixed, and the real-time filtered result is obtained using the newly measured and historical data.Figure 2
Schematic diagram of the sliding window.


## 3. Experimental Set-Up

### 3.1. Scaled Experimental Model

To validate the proposed method based on wavelet filtering, we used a 1/32 scaled model of an offshore supply vessel, as shown in Figure3, appended with a skeg and a bilge keel. The details of both the scaled model and the full-scale vessel are presented in Table 1. The vessel is equipped with three azimuth thrusters (two located aft and one in the fore) and one tunnel thruster (in the fore).Figure 3
1/32 scaled model of an offshore supply vessel.Table 1
Details of the offshore supply vessel.
ParameterSymbolUnitShipModelOverall lengthLOAm126.23.9438BreadthBm25.20.7875DepthDm10.60.3313Design draftTm7.00.2188DisplacementΔton13973.80.4264Experiments were conducted under winds, currents, and waves in the maneuvering and seakeeping basin of the China Ship Scientific Research Center (CSSRC). The basin is 69 m long and 46 m wide, and the water depth is 4 m. The current and wind load coefficients were obtained through tests in the wind tunnel. All the loads were simulated (the wind velocity was defined as 13 kn, and the current velocity was defined as 1.5 kn), and disturbances in terms of irregular waves with a significant wave height of 4.0 m and a characteristic period of 8.3 s were applied (corresponding to full-scale ones). The ITTC spectrum was employed in the model test to simulate the wave environment. Waves were generated by a wave generator with a system of 188 rocker flaps. Wind and current were generated by a large number of wind fans and a water pump system, respectively.The ITTC spectrum is defined as follows:(13)Sζ=173H1/32ω5T014exp−691ω4T014,where H1/3 is the significant wave height, T01 is the characteristic period, and ω is the circular wave frequency.The irregular seas were adjusted such that the spectral density distribution of the generated sea is compared well with the required theoretical energy distribution before the experiment. The requirements for the wave spectrum were set as follows: the significant wave height (H1/3) and the spectral characteristic period (T01) were within 5% of the theoretical values. The measured and target wave spectra are shown in Figure 4.Figure 4
Measured and target wave spectra.


### 3.2. Full-Scale Experimental Vessel

To verify the ability of wavelet filtering for the full-scale vessel, full-scale experiments were conducted on EXPLOERE I, a scientific investigation vessel, with the DP system (the DP control system was developed by CSSRC and approved by the China Classification Society), as shown in Figure5. The details of EXPLOERE I are presented in Table 2. The vessel was fitted with two main propellers and four tunnel thrusters (two at the stern and two in the bow).Figure 5
EXPLORER I scientific investigation vessel with DP system: (a) EXPLORER I and (b) CSSRC DP system.
(a)(b)Table 2
Details of the EXPLORER I scientific investigation vessel.
ParameterSymbolUnitShipOverall lengthLOAm89.6BreadthBm17.9DepthDm10.6Design draftTm5.5DisplacementΔton5533.0The CSSRC DP control system uses data from gyrocompasses, differential Global Positioning System (DGPS), mechanical reference systems (taut wire), and the BeiDou satellite position reference system, which allow the DP control system to position the vessel.The operating condition was sea state 4, the significant wave height and main wave period were measured to be 2.0 m and 6.0 s, respectively, by a wave radar, and the mean wind speed was measured to be 8–11 m/s by a wind sensor. The surface current speed was not measured in this experiment.


## 3.1. Scaled Experimental Model

To validate the proposed method based on wavelet filtering, we used a 1/32 scaled model of an offshore supply vessel, as shown in Figure3, appended with a skeg and a bilge keel. The details of both the scaled model and the full-scale vessel are presented in Table 1. The vessel is equipped with three azimuth thrusters (two located aft and one in the fore) and one tunnel thruster (in the fore).Figure 3
1/32 scaled model of an offshore supply vessel.Table 1
Details of the offshore supply vessel.
ParameterSymbolUnitShipModelOverall lengthLOAm126.23.9438BreadthBm25.20.7875DepthDm10.60.3313Design draftTm7.00.2188DisplacementΔton13973.80.4264Experiments were conducted under winds, currents, and waves in the maneuvering and seakeeping basin of the China Ship Scientific Research Center (CSSRC). The basin is 69 m long and 46 m wide, and the water depth is 4 m. The current and wind load coefficients were obtained through tests in the wind tunnel. All the loads were simulated (the wind velocity was defined as 13 kn, and the current velocity was defined as 1.5 kn), and disturbances in terms of irregular waves with a significant wave height of 4.0 m and a characteristic period of 8.3 s were applied (corresponding to full-scale ones). The ITTC spectrum was employed in the model test to simulate the wave environment. Waves were generated by a wave generator with a system of 188 rocker flaps. Wind and current were generated by a large number of wind fans and a water pump system, respectively.The ITTC spectrum is defined as follows:(13)Sζ=173H1/32ω5T014exp−691ω4T014,where H1/3 is the significant wave height, T01 is the characteristic period, and ω is the circular wave frequency.The irregular seas were adjusted such that the spectral density distribution of the generated sea is compared well with the required theoretical energy distribution before the experiment. The requirements for the wave spectrum were set as follows: the significant wave height (H1/3) and the spectral characteristic period (T01) were within 5% of the theoretical values. The measured and target wave spectra are shown in Figure 4.Figure 4
Measured and target wave spectra.


## 3.2. Full-Scale Experimental Vessel

To verify the ability of wavelet filtering for the full-scale vessel, full-scale experiments were conducted on EXPLOERE I, a scientific investigation vessel, with the DP system (the DP control system was developed by CSSRC and approved by the China Classification Society), as shown in Figure5. The details of EXPLOERE I are presented in Table 2. The vessel was fitted with two main propellers and four tunnel thrusters (two at the stern and two in the bow).Figure 5
EXPLORER I scientific investigation vessel with DP system: (a) EXPLORER I and (b) CSSRC DP system.
(a)(b)Table 2
Details of the EXPLORER I scientific investigation vessel.
ParameterSymbolUnitShipOverall lengthLOAm89.6BreadthBm17.9DepthDm10.6Design draftTm5.5DisplacementΔton5533.0The CSSRC DP control system uses data from gyrocompasses, differential Global Positioning System (DGPS), mechanical reference systems (taut wire), and the BeiDou satellite position reference system, which allow the DP control system to position the vessel.The operating condition was sea state 4, the significant wave height and main wave period were measured to be 2.0 m and 6.0 s, respectively, by a wave radar, and the mean wind speed was measured to be 8–11 m/s by a wind sensor. The surface current speed was not measured in this experiment.


## 4. Results and Discussion

### 4.1. Scaled-Model Experiments

Prior to the experimental verification of wavelet filtering in DP, experiments on a scaled model were conducted with EKF. The three degrees of freedom of the horizontal plan (sway, surge, and heading) were analyzed corresponding to the wave direction of 180° (coming from the vessel bow) by using wavelet filtering. The motion of the scaled model was measured using an optical measurement system.To test the effectiveness of the proposed method, we employed the stand deviation (SD)σx and the signal-to-noise ratio (SNR) in order to evaluate the reconstruction accuracy of real signals:(14)σx=1/N∑i=1Nxi−x¯2,SNR=10log∑i=1Nxi2/xi−x^i2,where xi denotes the measured signal and x^i denotes the filtered signal.The wavelet function, sliding window length, decomposition level, and thresholding method directly affect the filtering performance. In this study, we tested a range of window lengths and near-symmetric symlet wavelets (sym2, sym4, sym6, and sym8), as shown in Table3. The optimal decomposition level was 4, and the detailed signals of levels 1–3 were processed using the hard thresholding method. The SD and SNR do not change when the window size is ≥26 for different symlet wavelets. To validate the effectiveness of wavelet filtering, EKF has been used in [7]. The SD and SNR of EKF are 0.0335 and 8.6018, respectively, which is closer to the filtering performance of the symN wavelet (N = 4, 6, and 8), as shown in Table 3. Considering that higher-order N of vanishing moments will make the filtering signal smoother while consuming more time for filtering, the sym6 wavelet was selected for wave filtering in this model experiment.Table 3
SD and SNR of surge motion for different wavelets and sliding window lengths.
SymNWindow lengthWavelet filteringEKFSDSNRSDSNRSym2240.027710.29330.03358.6018250.027810.2631260.027810.2646270.027810.2646280.027810.2646290.027810.26462100.027810.2646Sym4240.03538.1797250.03338.6897260.03308.7648270.03308.7647280.03308.7647290.03308.76472100.03308.7647Sym6240.03608.5351250.03398.5929260.03378.5968270.03378.5968280.03378.5968290.03378.59682100.03378.5968Sym8240.03647.9172250.03398.5442260.03378.6012270.03368.6083280.03368.6083290.03368.60832100.03368.6083Furthermore, the relevant source code for wave filtering was developed using the LabVIEW platform on a personal computer with an Intel Core i7-4790 processor. To validate the program for real-time filtering, the time consumed can be determined using the LabVIEW Tick Count function. By filtering 1500 samples, different sliding window lengths with the sym6 wavelet were detected, and the average and maximum time consumed were obtained, as shown in Table4. It can be concluded that when the window length is 26, the execution time is around 0.4 ms, and the time consumed is around 15 times greater than that when the window length is 210. The sample time is 200 ms in the model experiment; thus, real-time filtering can be realized using the sym6 wavelet.Table 4
Execution time of wave filtering based on sym6 wavelet with different sliding window lengths.
Wave directionSymlet waveletStatisticExecute time (ms)Window length26272829210180°sym6Average0.360.721.443.005.84Maximum0.410.781.533.356.08As mentioned above, the smaller the window size, the shorter is the time consumed. However, the effectiveness of wave filtering must also be guaranteed. From the results presented in Tables3 and 4, we comprehensively considered the computational complexity and filtering performance for different wavelet functions, and we selected the sym6 wavelet with a window length of 26 for appropriate wavelet filtering in the model experiments.Furthermore, the reasons for selecting the optimal decomposition level of 4 and processing the detailed signals of levels 1–3 by the hard thresholding method will be analyzed.As shown in Figure4, the WF components corresponding to the wave spectrum are nearly in the interval [0.5 3]. The high-pass and low-pass filters are complementary and have a cut-off frequency of fs/2 at each level in the wavelet optimal decomposition. When the level of decomposition is 4, the theoretical approximation coefficients of the fourth level are in [0 0.3125] (the sample frequency fs is 5 Hz in this DP experiment).The decomposed approximation signals and detailed signals for surge motion are separately given with level 1–4 decomposition using the sym6 wavelet, as shown in Figure6. Further, the decomposed signals are subjected to fast Fourier transformation (FFT) and amplitude (AM) spectrum analysis, as shown in Figure 7. As can be seen, the corresponding high cut-off frequency of the detailed signals is around 0.3 in the fourth level. Therefore, the motion of the WF components can filter out the detailed signals of levels 1–3 by the hard thresholding method. Figure 8 shows the measured and filtered results of surge motion; the time evolution of the WF and LF components of the motion and their estimates are presented by wavelet filtering with SBE, wavelet filtering without SBE, and EKF. As can be seen, the second-order LF motion induced by slowly varying disturbances is effectively separated from the first-order WF oscillatory motion using wavelet filtering with SBE. For a clearer comparison, Figure 8(b) shows a magnified view of Figure 8(a). As can be seen, wavelet filtering provides a smoother estimation of the LF part of the motion compared with that obtained using EKF; the end effects of the filtered signal are restrained with SBE compared with that without SBE; the phase lag and amplitude distortion are insignificant because of the near-symmetric sym6 wavelet; and the delay time is around 0.4 s.Figure 6
Four decomposed levels of surge motion in the wave direction of 180° with the sym6 wavelet: (a) levels 1–4 corresponding to approximation signal xa and (b) levels 1–4 corresponding to detailed signal xd.
(a)(b)Figure 7
Amplitude (AM) spectrum of surge motion in the wave direction of 180° with levels 1–4: (a) FFT of approximation signal xa and (b) FFT of detailed signal xd.
(a)(b)Figure 8
Experimental results of surge motion in the wave direction of 180°. Total motion and estimation of the LF components of motion with the wavelet transform and EKF: (a) surge motion and (b) magnified view of (a). The body-frame xyz is fixed to the model along the path of the model, with the positivex-axis (surge) in the forward direction and the positive y-axis (sway) to the starboard.
(a)(b)To further verify the effectiveness of real-time wave filtering with the sym6 wavelet in different wave directions, DP station keeping experiments were conducted in the wave direction of 135° on the basis of the aforementioned selections of the wavelet function, window length, decomposition level, and hard thresholding method. The sway, surge motion, and heading were analyzed by real-time wavelet filtering, as shown in Figures9–11. As can be seen, the same effective filtering performance is verified compared to EKF.Figure 9
Experimental results of surge motion in the wave direction of 135°. Measured motion and estimation of the LF components of motion with the wavelet transform and EKF: (a) sway motion and (b) the magnified view of (a).
(a)(b)Figure 10
Experimental results of sway motion in the wave direction of 135°. Measured motion and estimation of the LF components of motion with the wavelet transform and EKF: (a) sway motion and (b) the magnified view of (a).
(a)(b)Figure 11
Experimental results of heading in the wave direction of 135°. Measured motion and estimation of the LF components of motion with the wavelet transform and EKF: (a) heading motion and (b) the magnified view of (a).
(a)(b)


### 4.2. Full-Scale Experiment

To verify the feasibility of wavelet filtering in a full-scale vessel, first, the offline data were analyzed by conducting a full-scale experiment on the EXPLOERE I scientific investigation vessel on the basis of the wavelet filtering with SBE and the sliding window method. Compared to the scaled-model experiment, the wave spectral density distribution was not accurately measured in the full-scale experiment. Hence, the decomposition level and thresholding method were selected by amplitude spectrum analysis of the offline data.The measured signal frequency was distributed on the whole frequency coordinate axis, and the amplitude spectrum distribution was in the high-frequency components, as shown in Figure12. Clearly, the amplitude spectrum distribution with frequencies greater than 0.5 Hz was reduced to 0. Through a comparison of the amplitude spectra of the measured signal and the filtering signal, we found that the signal energy after filtering was clearly reduced to 0 in the high-frequency components. Hence, the decomposition level of 4 and the detailed signals of levels 1–3 were processed by the hard thresholding method and selected in the full-scale experiments.Figure 12
Amplitude spectra of measured and offline wavelet filtering motion for north, east and heading. The earth-fixed reference XEYEZE is denoted with the positive XE-axis to the north and the positive y-axis to the east.The same sym6 wavelet and the sliding window length of 26 were used as a reference for real-time filtering in the full-scale experiments. Then, the high- and low-frequency components after threshold filtering were combined for signal reconstruction. The final filtering signal was obtained through reconstruction, as shown in Figure 13. Thus, the LF signal waveform was effectively recovered after wavelet filtering.Figure 13
Experimental results of north, east and heading: measured motion and offline estimation of the LF components of motion with the wavelet filter.According to the aforementioned results, we could further verify the effectiveness of wavelet filtering. Online data filtering was conducted in the full-scale experiments on the basis of real-time wavelet filtering with SBE and the sliding window method. The experimental results are shown in Figure14. For the real-time filtering results, the amplitude spectrum analysis is shown in Figure 15. As can be seen, the WF motion component was effectively filtered out and the smoother LF motion waveform was effectively recovered after filtering.Figure 14
Experimental results of north, east and heading: measured motion and online real-time estimation of the LF components of motion with the wavelet filter, SBE and sliding window method.Figure 15
Amplitude spectra of measured and online wavelet filtered motion for north, east and heading.


### 4.3. Comparison of Filtering Performance

Furthermore, we quantitatively compared the filtering performance of wavelet filtering and EKF. The filtering approaches were compared and the cost functionE is used as follows:(15)Eωc=∫ωcωs/2sWFωdω∫0ωs/2smωdω×100%,where ωc is a cut-off frequency of WF component, ωs is circular frequency corresponding to sample frequency fs, sm(ω) is energy spectral density function of the measured signal, and sWF(ω) is energy spectral density function of the measured signal of WF motion. E represents the energy with the frequency band ωc≤ω≤ωs; the smaller the value E, the better the filtering performance.The energy spectral density functions of the model and full-scale experimental results are shown in Figures16–18. According to the results, we can see that the low-frequency motion components are effectively separated, and wave-frequency motion components can be filtered out by wavelet filtering.Figure 16
Energy spectral density function of the model experimental results in wave direction 180°: (a) surge (x direction) motion, (b) sway (y direction) motion and (c) heading.
(a)(b)(c)Figure 17
Energy spectral density function of the model experimental results in wave direction 135°: (a) surge (x direction) motion; (b) sway (y direction) motion and (c) heading.
(a)(b)(c)Figure 18
Energy spectral density function of the full-scale experimental results: (a) north (x direction), (b) east (y direction) and (c) heading.
(a)(b)(c)The filtering methods are compared usingE(ωc) for the aforementioned scaled model and full-scale vessel, and three horizontal degrees of freedom (X or surge, Y or sway, and heading) were analyzed (ωc is 0.625 π·rad/s), and horizontal three degrees of freedom motion (x or north, y or east, and heading) are analyzed. The results of E are shown in Table 5. As can be seen from Table 5, compared to the EKF approach, the energy of wave-frequency components is lower using wavelet filtering, and the wavelet filtering performs better for all cases in the model and full-scale experiments.Table 5
Performance comparison between wavelet filtering and EKF.
Wave directionMotionsEWF-wavelet (%)EWF-EKF (%)EWF-m (%)Scaled180°x0.4071.30012.269y0.0210.0750.176Heading0.0890.3003.954135°x0.1110.6524.489y0.0040.0310.222Heading0.0440.2712.392Full-scale—x0.3180.6114.953y0.4470.9585.210Heading0.2170.3141.658As mentioned above, compared with the EKF with ship dynamics, the wavelet filtering returns better results of position and heading filtering. However, the most important point we presented in the work is that the wavelet filtering approach is model-independent, and the low-frequency and wave-frequency motion model is not established compared to the EKF, especially a damped oscillation wave-frequency model is not accurate in EKF, and the characteristics of waves still need to be considered to obtain the better filtering performance.EWF-wavelet is the WF motion using wavelet filtering; EWF-EKF is the WF motion using EKF; and EWF-m is the WF components of measured motion.


## 4.1. Scaled-Model Experiments

Prior to the experimental verification of wavelet filtering in DP, experiments on a scaled model were conducted with EKF. The three degrees of freedom of the horizontal plan (sway, surge, and heading) were analyzed corresponding to the wave direction of 180° (coming from the vessel bow) by using wavelet filtering. The motion of the scaled model was measured using an optical measurement system.To test the effectiveness of the proposed method, we employed the stand deviation (SD)σx and the signal-to-noise ratio (SNR) in order to evaluate the reconstruction accuracy of real signals:(14)σx=1/N∑i=1Nxi−x¯2,SNR=10log∑i=1Nxi2/xi−x^i2,where xi denotes the measured signal and x^i denotes the filtered signal.The wavelet function, sliding window length, decomposition level, and thresholding method directly affect the filtering performance. In this study, we tested a range of window lengths and near-symmetric symlet wavelets (sym2, sym4, sym6, and sym8), as shown in Table3. The optimal decomposition level was 4, and the detailed signals of levels 1–3 were processed using the hard thresholding method. The SD and SNR do not change when the window size is ≥26 for different symlet wavelets. To validate the effectiveness of wavelet filtering, EKF has been used in [7]. The SD and SNR of EKF are 0.0335 and 8.6018, respectively, which is closer to the filtering performance of the symN wavelet (N = 4, 6, and 8), as shown in Table 3. Considering that higher-order N of vanishing moments will make the filtering signal smoother while consuming more time for filtering, the sym6 wavelet was selected for wave filtering in this model experiment.Table 3
SD and SNR of surge motion for different wavelets and sliding window lengths.
SymNWindow lengthWavelet filteringEKFSDSNRSDSNRSym2240.027710.29330.03358.6018250.027810.2631260.027810.2646270.027810.2646280.027810.2646290.027810.26462100.027810.2646Sym4240.03538.1797250.03338.6897260.03308.7648270.03308.7647280.03308.7647290.03308.76472100.03308.7647Sym6240.03608.5351250.03398.5929260.03378.5968270.03378.5968280.03378.5968290.03378.59682100.03378.5968Sym8240.03647.9172250.03398.5442260.03378.6012270.03368.6083280.03368.6083290.03368.60832100.03368.6083Furthermore, the relevant source code for wave filtering was developed using the LabVIEW platform on a personal computer with an Intel Core i7-4790 processor. To validate the program for real-time filtering, the time consumed can be determined using the LabVIEW Tick Count function. By filtering 1500 samples, different sliding window lengths with the sym6 wavelet were detected, and the average and maximum time consumed were obtained, as shown in Table4. It can be concluded that when the window length is 26, the execution time is around 0.4 ms, and the time consumed is around 15 times greater than that when the window length is 210. The sample time is 200 ms in the model experiment; thus, real-time filtering can be realized using the sym6 wavelet.Table 4
Execution time of wave filtering based on sym6 wavelet with different sliding window lengths.
Wave directionSymlet waveletStatisticExecute time (ms)Window length26272829210180°sym6Average0.360.721.443.005.84Maximum0.410.781.533.356.08As mentioned above, the smaller the window size, the shorter is the time consumed. However, the effectiveness of wave filtering must also be guaranteed. From the results presented in Tables3 and 4, we comprehensively considered the computational complexity and filtering performance for different wavelet functions, and we selected the sym6 wavelet with a window length of 26 for appropriate wavelet filtering in the model experiments.Furthermore, the reasons for selecting the optimal decomposition level of 4 and processing the detailed signals of levels 1–3 by the hard thresholding method will be analyzed.As shown in Figure4, the WF components corresponding to the wave spectrum are nearly in the interval [0.5 3]. The high-pass and low-pass filters are complementary and have a cut-off frequency of fs/2 at each level in the wavelet optimal decomposition. When the level of decomposition is 4, the theoretical approximation coefficients of the fourth level are in [0 0.3125] (the sample frequency fs is 5 Hz in this DP experiment).The decomposed approximation signals and detailed signals for surge motion are separately given with level 1–4 decomposition using the sym6 wavelet, as shown in Figure6. Further, the decomposed signals are subjected to fast Fourier transformation (FFT) and amplitude (AM) spectrum analysis, as shown in Figure 7. As can be seen, the corresponding high cut-off frequency of the detailed signals is around 0.3 in the fourth level. Therefore, the motion of the WF components can filter out the detailed signals of levels 1–3 by the hard thresholding method. Figure 8 shows the measured and filtered results of surge motion; the time evolution of the WF and LF components of the motion and their estimates are presented by wavelet filtering with SBE, wavelet filtering without SBE, and EKF. As can be seen, the second-order LF motion induced by slowly varying disturbances is effectively separated from the first-order WF oscillatory motion using wavelet filtering with SBE. For a clearer comparison, Figure 8(b) shows a magnified view of Figure 8(a). As can be seen, wavelet filtering provides a smoother estimation of the LF part of the motion compared with that obtained using EKF; the end effects of the filtered signal are restrained with SBE compared with that without SBE; the phase lag and amplitude distortion are insignificant because of the near-symmetric sym6 wavelet; and the delay time is around 0.4 s.Figure 6
Four decomposed levels of surge motion in the wave direction of 180° with the sym6 wavelet: (a) levels 1–4 corresponding to approximation signal xa and (b) levels 1–4 corresponding to detailed signal xd.
(a)(b)Figure 7
Amplitude (AM) spectrum of surge motion in the wave direction of 180° with levels 1–4: (a) FFT of approximation signal xa and (b) FFT of detailed signal xd.
(a)(b)Figure 8
Experimental results of surge motion in the wave direction of 180°. Total motion and estimation of the LF components of motion with the wavelet transform and EKF: (a) surge motion and (b) magnified view of (a). The body-frame xyz is fixed to the model along the path of the model, with the positivex-axis (surge) in the forward direction and the positive y-axis (sway) to the starboard.
(a)(b)To further verify the effectiveness of real-time wave filtering with the sym6 wavelet in different wave directions, DP station keeping experiments were conducted in the wave direction of 135° on the basis of the aforementioned selections of the wavelet function, window length, decomposition level, and hard thresholding method. The sway, surge motion, and heading were analyzed by real-time wavelet filtering, as shown in Figures9–11. As can be seen, the same effective filtering performance is verified compared to EKF.Figure 9
Experimental results of surge motion in the wave direction of 135°. Measured motion and estimation of the LF components of motion with the wavelet transform and EKF: (a) sway motion and (b) the magnified view of (a).
(a)(b)Figure 10
Experimental results of sway motion in the wave direction of 135°. Measured motion and estimation of the LF components of motion with the wavelet transform and EKF: (a) sway motion and (b) the magnified view of (a).
(a)(b)Figure 11
Experimental results of heading in the wave direction of 135°. Measured motion and estimation of the LF components of motion with the wavelet transform and EKF: (a) heading motion and (b) the magnified view of (a).
(a)(b)


## 4.2. Full-Scale Experiment

To verify the feasibility of wavelet filtering in a full-scale vessel, first, the offline data were analyzed by conducting a full-scale experiment on the EXPLOERE I scientific investigation vessel on the basis of the wavelet filtering with SBE and the sliding window method. Compared to the scaled-model experiment, the wave spectral density distribution was not accurately measured in the full-scale experiment. Hence, the decomposition level and thresholding method were selected by amplitude spectrum analysis of the offline data.The measured signal frequency was distributed on the whole frequency coordinate axis, and the amplitude spectrum distribution was in the high-frequency components, as shown in Figure12. Clearly, the amplitude spectrum distribution with frequencies greater than 0.5 Hz was reduced to 0. Through a comparison of the amplitude spectra of the measured signal and the filtering signal, we found that the signal energy after filtering was clearly reduced to 0 in the high-frequency components. Hence, the decomposition level of 4 and the detailed signals of levels 1–3 were processed by the hard thresholding method and selected in the full-scale experiments.Figure 12
Amplitude spectra of measured and offline wavelet filtering motion for north, east and heading. The earth-fixed reference XEYEZE is denoted with the positive XE-axis to the north and the positive y-axis to the east.The same sym6 wavelet and the sliding window length of 26 were used as a reference for real-time filtering in the full-scale experiments. Then, the high- and low-frequency components after threshold filtering were combined for signal reconstruction. The final filtering signal was obtained through reconstruction, as shown in Figure 13. Thus, the LF signal waveform was effectively recovered after wavelet filtering.Figure 13
Experimental results of north, east and heading: measured motion and offline estimation of the LF components of motion with the wavelet filter.According to the aforementioned results, we could further verify the effectiveness of wavelet filtering. Online data filtering was conducted in the full-scale experiments on the basis of real-time wavelet filtering with SBE and the sliding window method. The experimental results are shown in Figure14. For the real-time filtering results, the amplitude spectrum analysis is shown in Figure 15. As can be seen, the WF motion component was effectively filtered out and the smoother LF motion waveform was effectively recovered after filtering.Figure 14
Experimental results of north, east and heading: measured motion and online real-time estimation of the LF components of motion with the wavelet filter, SBE and sliding window method.Figure 15
Amplitude spectra of measured and online wavelet filtered motion for north, east and heading.


## 4.3. Comparison of Filtering Performance

Furthermore, we quantitatively compared the filtering performance of wavelet filtering and EKF. The filtering approaches were compared and the cost functionE is used as follows:(15)Eωc=∫ωcωs/2sWFωdω∫0ωs/2smωdω×100%,where ωc is a cut-off frequency of WF component, ωs is circular frequency corresponding to sample frequency fs, sm(ω) is energy spectral density function of the measured signal, and sWF(ω) is energy spectral density function of the measured signal of WF motion. E represents the energy with the frequency band ωc≤ω≤ωs; the smaller the value E, the better the filtering performance.The energy spectral density functions of the model and full-scale experimental results are shown in Figures16–18. According to the results, we can see that the low-frequency motion components are effectively separated, and wave-frequency motion components can be filtered out by wavelet filtering.Figure 16
Energy spectral density function of the model experimental results in wave direction 180°: (a) surge (x direction) motion, (b) sway (y direction) motion and (c) heading.
(a)(b)(c)Figure 17
Energy spectral density function of the model experimental results in wave direction 135°: (a) surge (x direction) motion; (b) sway (y direction) motion and (c) heading.
(a)(b)(c)Figure 18
Energy spectral density function of the full-scale experimental results: (a) north (x direction), (b) east (y direction) and (c) heading.
(a)(b)(c)The filtering methods are compared usingE(ωc) for the aforementioned scaled model and full-scale vessel, and three horizontal degrees of freedom (X or surge, Y or sway, and heading) were analyzed (ωc is 0.625 π·rad/s), and horizontal three degrees of freedom motion (x or north, y or east, and heading) are analyzed. The results of E are shown in Table 5. As can be seen from Table 5, compared to the EKF approach, the energy of wave-frequency components is lower using wavelet filtering, and the wavelet filtering performs better for all cases in the model and full-scale experiments.Table 5
Performance comparison between wavelet filtering and EKF.
Wave directionMotionsEWF-wavelet (%)EWF-EKF (%)EWF-m (%)Scaled180°x0.4071.30012.269y0.0210.0750.176Heading0.0890.3003.954135°x0.1110.6524.489y0.0040.0310.222Heading0.0440.2712.392Full-scale—x0.3180.6114.953y0.4470.9585.210Heading0.2170.3141.658As mentioned above, compared with the EKF with ship dynamics, the wavelet filtering returns better results of position and heading filtering. However, the most important point we presented in the work is that the wavelet filtering approach is model-independent, and the low-frequency and wave-frequency motion model is not established compared to the EKF, especially a damped oscillation wave-frequency model is not accurate in EKF, and the characteristics of waves still need to be considered to obtain the better filtering performance.EWF-wavelet is the WF motion using wavelet filtering; EWF-EKF is the WF motion using EKF; and EWF-m is the WF components of measured motion.


## 5. Conclusions

A real-time wave filtering approach for DP systems was proposed on the basis of the wavelet transform, SBE, and the sliding window method. We analyzed the performance of wavelet filtering using DWT for separating the high-frequency oscillatory motion from the second-order low-frequency slowly varying motion by using a scaled model of an offshore supply vessel and conducting full-scale experiments on the EXPLOERE I scientific investigation vessel. Specifically, the noise component was in the high-frequency part of the wavelet coefficient, the hard thresholding method was preferred, and a decomposition level of 4 was used for wavelet filtering. The scaled model and full-scale experimental results demonstrated that the proposed approach can filter out the high frequencies of a raw signal by retaining only the low frequencies, and no significant signal delay is caused by the selection of a near-symmetric wavelet function. In addition, different sliding window sizes were analyzed, and the end effects with and without SBE were compared. By comprehensively considering the computational complexity and filtering performance for different wavelet families, we selected the sym6 wavelet and a window length of 26 for appropriate wavelet filtering in the scaled model and full-scale experiments. The proposed real-time wavelet filtering method was shown to be suitable for nonstationary and nonlinear signals of the position and heading of DP vessels, containing high-frequency components. Moreover, it was found to significantly improve the quality of the separated components without relying on the previous knowledge of dynamics of the vessel.Although the current performance of real-time wavelet filtering seems acceptable in DP experiments on the scaled model and full-scale vessel, some problems remain owing to the different characteristics of the wave spectrum in different sea operation areas for a full-scale vessel. These problems must be addressed in future studies. Thus, further investigation is required to fully develop the adaptive wavelet filtering approach in different sea states. Moreover, considering that the motion periods will be different for different types of vessels under waves, to realize a wide range of engineering applications for DP vessels, the wavelet filtering approach should be further validated using different types of vessels.



---

*Source: 5445420-2022-11-22.xml*