A new controller design of electrohydraulic servo system based on empirical mode decomposition
Jing Huang^{1} , Changchun Li^{2} , Hao Yan^{3} , Xuesong Yang^{4} , Jing Li^{5}
^{1, 2, 3, 4, 5}Beijing Jiaotong University, Beijing, China
^{1}Corresponding author
Journal of Vibroengineering, Vol. 18, Issue 2, 2016, p. 927937.
Received 9 November 2015; received in revised form 25 February 2016; accepted 23 March 2016; published 31 March 2016
JVE Conferences
The signal of electrohydraulic servo system is nonstationary and timevarying due to the influence of vibration, noise and mechanical impact. The traditional digital filter always suffers delay in time domain and the delay increases along with the increasing of frequency. Considering the features of electrohydraulic servo system, the HilbertHuang transform method is an effective method to decompose the original signal and obtain the noise components. Some improvements are made based on Hilbert Huang transform method and a new real time online filtering method is proposed in this paper. This improved filter is able to decompose out the noise components and other interference components from original signal, and remove them off in real time. Based on this new online filter, a new controller is also designed. Compared the filtering result with the traditional digital filter, this new controller’s control performance is much better.
Keywords: online filter, HHT, real time, electrohydraulic servo system.
1. Introduction
Due to the influence of the hydraulic transmission medium’s flow ability, compressibility, viscosity, and the influence of working environment such as temperature and pressure, the hydraulic cylinder sometimes works at some unusual states such as vibration, noise and impact [1]. The signal of electrohydraulic servo system is nonstationary and timevariation. To filter these noise signals, the general solution is adding a digital filter into system. However, this solution also brings some problems such as phase delay and amplitude attenuation.
HilbertHuang transform (HHT) method [2] is proposed by Norden E. Huang et al. in 1998. This method is suitable for nonlinear and nonstationary signal analysis. At present, HHT technology has been used widely in the field of fault diagnosis, signal processing [35], oceanography [6] and so on [7, 8], and has obtained some good results. Considering the characters of electrohydraulic servo system, HHT method can decompose out the noise components and other interference components from original signal, therefore, we are able to remove any unnecessary noise component. But the empirical mode decomposition (EMD) method which is a part of HHT method has some disadvantages when using in real time online filter, such as end point effect and computation efficiency. If we want to use the HHT method in electro hydraulic servo control system, some improvements on the EMD method are essential. Also, it is necessary to design a new controller based on this improved method.
2. Improvements of EMD method
2.1. Brief introduction of EMD process
EMD method can decompose out every single intrinsic mode function (IMF) from high frequency to low frequency, which is selfadaption. Most of these IMFs have their own physical meaning. The process of EMD method is described in the following steps 13.
Step 1: Find out all the local maximum and minimum points. The local maximum and minimum points of original signal $x$ is defined in Eq. (1) and Eq. (2):
Step 2: Select a kind of interpolation method and employ it between all extreme points, then get the envelope curve ${x}_{max}\left(t\right)$ and ${x}_{min}\left(t\right)$.
Step 3: Calculate the difference $h\left(t\right)$ according to Eq. (3). If $h\left(t\right)$ does not satisfy Eq. (4), then repeat steps 13:
The first $h\left(t\right)$ satisfied Eq. (4) is the first IMF, let ${I}_{1}\left(t\right)=h\left(t\right)$, $r\left(t\right)={x\left(t\right)I}_{1}\left(t\right)$, $r\left(t\right)$ be the original signal of next process round. Repeat steps 13, the rest of other IMF${I}_{2}\left(t\right)$, ${I}_{3}\left(t\right)$, ${I}_{4}\left(t\right)$, …, ${I}_{n}\left(t\right)$ will be decomposed out. When $r\left(t\right)$ is a monotonic or constant series, this decomposition process will be terminated.
2.2. Optimizations of decomposition process
Review the whole process of EMD, the envelope curve is the key of decomposition. But the end point effect has significant influence on the decomposition [3], especially when the length of data is short. Many researchers have done a lot of work [912] to restrain the end point effect. To get a better filtering performance of short data, we make optimization of the EMD from four aspects as stated below.
Optimization 1: If the end point is the maximum or minimum point of the whole data section, then the end point will be treated as a local maximum or minimum point:
Optimization 2: The decomposition process may miss some turning points if the process follows the Eq. (1) and Eq. (2) strictly. And these missed turning points will reduce the effect of envelop curve as shown in Fig. 1. In this figure, ${P}_{2}$ and ${P}_{19}$ are two turning points, and they are not treated as local extreme points, so the envelop curves at point ${P}_{4}$ and ${P}_{18}$ are not complete: the max envelope curve is missed at point ${P}_{4}$ and the min envelop curve is missed at point ${P}_{18}$. At next process round, the value of point ${P}_{4}$ will change to ${P}_{4}^{\text{'}}$:
and the same thing will happen to point ${P}_{18}$.
The decomposition result of Fig. 1 is shown in Fig. 2. We notice that the noise at point ${P}_{4}$ and ${P}_{18}$ is missed. Along with the further development of this decomposition process, ${P}_{5}$ and ${P}_{16}$ will also become 0. This effect will diffuse to the data center, then the data will be changed from $\left[{x}_{1},{x}_{2},{x}_{3},\dots {,x}_{n2},{x}_{n1},{x}_{n}\right]$ to $\left[\mathrm{0,0},0,\dots 0,{x}_{n/2p},\dots ,{x}_{n/21},{x}_{n/2},{x}_{n/2+1},\dots {x}_{n/2+q},0,\dots \mathrm{0,0},0\right]$. This means that the filtering effect only works in the data’s center area. The decomposition result under this situation is shown in Fig. 3, from this figure it’s easy to see that all the noise components (IMF1IMF3) are not ideal.
Therefore, the turning points should be included in the definition of local extremely points. The maximum point of turning point is redefined as Eqs. (5)(6), and the minimum point of turning point is redefined as Eqs. (7)(8);
Fig. 1. Effect of turning point on envelope (a)
Fig. 2. Effect of turning point on envelope (b)
Fig. 3. The decomposition result without considering turning points
After extending the definitions of local maximum and minimum points, the new max points, min points and envelope curves in Fig. 1 are shown in Fig. 4 and Fig. 5. The envelope curves can be obtained by using three spline interpolation method. With considering the point ${P}_{2}$ and ${P}_{19}$, the max envelop curve and min envelope curve become complete. The value of point ${P}_{4}$ changes to ${P}_{4}^{\text{'}\text{'}}$: ${{P}_{4}^{\text{'}\text{'}}=P}_{4}{(P}_{4max}+{P}_{4min})/2\ne 0$. At next process round, the decomposition result is shown is Fig. 5. In Figs. 4 and 5, we can notice that the diffusion effect of data changing into 0 no longer exists.
After using these revised definitions, the decomposition result of the signal in Fig. 3 is shown in Fig. 6. In this result, IMF1 and IMF2 are the noise components, IMF3 is the useful signal. Compared Fig. 6 with Fig. 3, the new filter works better, which means the optimization 2 is necessary and effective.
Optimization 3: The envelop curve may have overshoot or undershoot when the span between two local extreme points is too larger than others. This situation is shown in Fig. 7. As shown in this figure, the span between point ${P}_{2}$ and ${P}_{53}$ is ${d}_{1}$, ${d}_{1}=\mathrm{}$53$$2$=$51. The other spans are ${d}_{2}$ to ${d}_{4}$, ${d}_{2}=\text{68}\text{53}=\text{15}$, ${d}_{3}=\text{77}\text{68}=\text{9}$,${d}_{4}=\text{97}\text{77}=\text{20}$, the average value of these spans is ${d}_{mean}$, ${d}_{mean}={(d}_{1}+{d}_{2}+{d}_{3}+{d}_{4})/4=$23.75, ${d}_{1}>{\text{1.5}d}_{mean}$. The undershoot (between point ${P}_{2}$ and ${P}_{53}$) of envelope curve will have a significant impact on the results. In this example, the value of local max points’ envelope curve is smaller than local min points’ envelope curve, which will introduce errors into final result. In this undershoot section, the envelope curve is ought to be the signal itself, the correct calculation process should be:
implying that there is no noise in the data of this section. But due to this undershoot, the calculation result becomes ${P}_{i}^{\text{'}}={P}_{i}({P}_{imax}+{P}_{imin})/2\ne 0$, indicating that the process will introduce some error noises into result.
Fig. 4. The envelope effect with considering turning points (a)
Fig. 5. The envelope effect with considering turning points (b)
Fig. 6. The decomposition result with considering turning points
In order to avoid this situation, interpolation method could be an outstanding solution. Firstly, the average span ${d}_{mean}$ of $[{d}_{1},{d}_{2},{d}_{3},\dots ,{d}_{m}]$ can be calculated by Eq. (9), in which ${d}_{i}$ represents the span between two local extreme points ${P}_{p}$ and ${P}_{q}$, ${d}_{i}=pq$. Secondly, find out the maximum span according to Eq. (10) and (11). Thirdly, insert points in the data section which has the maximum span. If the data is $\left[{x}_{1},{x}_{2},{x}_{3},\dots ,{x}_{N}\right]$ and local extreme points are $[{e}_{1},{e}_{2},{e}_{3},\dots ,{e}_{n+1}]$, and the max span is between ${e}_{i}$ and ${e}_{j}$, $1\le i<j\le n+1$, so the exteme points will change to $[{e}_{1},{e}_{2},{e}_{3},\dots ,{e}_{i},{{x}_{i+1\bullet {d}_{mean}},{x}_{i+2\bullet {d}_{mean}},\dots ,{x}_{i+t\bullet {d}_{mean}},{e}_{j},\dots ,e}_{n+1}]$ after interpolation, in which $t=ji/{d}_{mean}$, $t\in N$:
in which $N$ is the length of data, $n$ is the number of span.
After applying this interpolation method, the envelop curve in Fig. 7 has been improved which is shown in Fig. 8.
Fig. 7. The decomposition result after considering turning points
Fig. 8. The envelop curve after interpolation
Optimization 4: In order to improve calculation efficiency and reduce operation time, a termination condition Eq. (12) is added. Because in some cases, the noise components are just in the first few IMFs, the filtering result will be inaccurate if the decomposition process runs to the final step ignoring Eq. (12). For example, if the noise components are $[{IMF}_{1},\dots ,{IMF}_{i}]$, the other components $[{IMF}_{i+1},\dots ,{IMF}_{n}]$ are useful signals. The original signal is $s$:
If the filter takes $[{IMF}_{i+1},\dots ,{IMF}_{n}]$ as noises, then the signal after filtering will be ${s}^{\text{'}}$:
This situation is shown in Fig. 9, only $IM{F}_{1}$ is the noise component, $IM{F}_{2}$ is useful signal, and the filter takes the other IMFs as noise components, the filtering result will be all zeros.
To avoid this situation, Eq. (12) is employed to judge the similarity between the original signal and the decomposition component. If the decomposition result is similar to the original signal, the decomposition process should be stopped:
in which ${X}_{i}$ is the original signal, ${I}_{i}$ is the $IMF$ component.
Fig. 9. The noise component is only in IMF1
3. The Design of online filter
The filtering process is carried out according to the following steps 18.
Step 1: Set the length of filtering window as n, when system acquires a new data point ${x}_{n+1}$, update the data of this window from $\left[{x}_{1},{x}_{2},{x}_{3},\dots ,{x}_{n}\right]$ to $\left[{x}_{2},{x}_{3},{x}_{4},\dots ,{x}_{n+1}\right]$.
Step 2: Find out all local maximum points ${x}_{max}=\left[{x}_{max1},{x}_{max2},{x}_{max3},\dots ,{x}_{maxp}\right]$ according to the Eqs. (1), (5), (6), and all local minimum points ${x}_{min}=\left[{x}_{min1},{x}_{min2},{x}_{min3},\dots ,{x}_{minq}\right]$ according to the Eqs. (2), (7), (8). This step also follows the optimization 1.
Step 3: Following optimization 3, ${x}_{max}$ and ${x}_{min}$ become ${x}_{max}^{\text{'}}$and ${x}_{min}^{\text{'}}$.
Step 4: In the array of ${x}_{max}^{\text{'}}$ and ${x}_{min}^{\text{'}}$, the three spine interpolation method is performed in range $\left[max1,maxp\right]$ and $\left[min1,minq\right]$, so the envelop of maximum points is ${e}_{max}$:
In the same way, ${e}_{min}=\left[{e}_{min1},{e}_{min1+1},{e}_{min1+2},\dots ,{e}_{minq}\right]$. The rest sections in range $\left[1,n\right]$ take the original data itself as the envelope curve, after this the ${e}_{max}$ and ${e}_{min}$ become ${e}_{max}^{\text{'}}$ and ${e}_{min}^{\text{'}}$:
Step 5: Let $m\left(t\right)=\left({e}_{max}^{\text{'}}\left(t\right)+{e}_{min}^{\text{'}}\left(t\right)\right)/2$, then $h\left(t\right)=x\left(t\right)m\left(t\right)$. Determine whether the process needs to be terminated according to Eqs. (4) and (12). If not, repeat the steps 25 until Eqs. (4) and (12) are satisfied. After these steps, the rest of $h\left(t\right)$ is one $IMF$, i.e., $IM{F}_{1}$.
Step 6: The original data $x\left(t\right)$ minus ${IMF}_{1}$ is the original signal $r\left(t\right)$ of next process round. For signal $r\left(t\right)$, repeat steps 25 then ${IMF}_{2}$ can be decomposed out. In the same way, the others ${IMF}_{3}$, ${IMF}_{4}$,…, ${IMF}_{n}$ also can be decomposed out.
This whole decomposition process can be represented in Table 1.
3.1. The order of noise components
Before the start of decomposition process, the order of noise components needs to be determined. Decompose the acquired data of electrohydraulic servo system by using EMD method, and the result is shown in Fig. 10. The noise components are mainly in the first few components $IM{F}_{1}$$IM{F}_{5}$, therefore, the process should be ended at $IM{F}_{5}$. The signal $x\left(t\right)$ becomes $y\left(t\right)$ after filtering, $y\left(t\right)=x\left(t\right){IMF}_{1}{IMF}_{2}{IMF}_{3}{IMF}_{4}{IMF}_{5}$.
Table 1. The whole process of filtering decomposition
Signal

$IM{F}_{1}$

$IM{F}_{2}$

$IM{F}_{3}$

$\cdots $

$IM{F}_{n}$

$X$

$x{IMF}_{1}={r}_{1}$

${r}_{1}{IMF}_{2}={r}_{2}$

${r}_{n2}{IMF}_{n1}={r}_{n1}$


0

$x{m}_{1}={h}_{1}$

${r}_{1}$

${m}_{2}={h}_{2}$

${r}_{2}$

${m}_{3}={h}_{3}$

${r}_{n1}$

${m}_{n}={h}_{n}$


1

${h}_{1}$

${m}_{11}={h}_{11}$

${h}_{2}$

${m}_{21}={h}_{21}$

${h}_{3}$

${m}_{31}={h}_{31}$

${h}_{n}$

${m}_{n1}={h}_{n1}$


2

${h}_{11}$

${m}_{12}={h}_{12}$

${h}_{21}$

${m}_{22}={h}_{22}$

${h}_{31}$

${m}_{32}={h}_{32}$

${h}_{n1}$

${m}_{n2}={h}_{n2}$


3

${h}_{12}$

${m}_{13}={h}_{13}$

${h}_{22}$

${m}_{23}={h}_{23}$

${h}_{32}$

${m}_{33}={h}_{33}$

${h}_{n2}$

${m}_{n3}={h}_{n3}$


$\vdots $

$\vdots $

$\vdots $

$\vdots $

$\vdots $


$k$

${h}_{1(k1)}{m}_{1k}={h}_{1k}$

${h}_{2(k1)}{m}_{2k}={h}_{2k}$

${h}_{3(k1)}{m}_{3k}={h}_{3k}$

${h}_{n(k1)}{m}_{nk}={h}_{nk}$


$k+$1

${IMF}_{1}={h}_{1k}$

${IMF}_{2}={h}_{2k}$

${IMF}_{3}={h}_{3k}$

${IMF}_{n}={h}_{nk}$

Fig. 10. The decomposition result
3.2. The length of filtering window
The length of filtering window is vital to filtering process. On one side, the longer the window is, the more information it contains. On the other side, the longer the window is, the more time and the bigger phase delay it takes. Therefore, a tradeoff between these two sides is necessary. Compare the results of length 20 and length 50, which are shown in Fig. 11 and Fig. 12, respectively, the length of this system is determined to 50.
3.3. The properties of the filter
Property 1: This new filter is selfadaptive, because this filter works in time domain and the decomposition process is based on the characteristic of the signal itself. When we design this new filter, the interpolation method, the orders of noise and the length of filter window are priorities need to be considered, the frequency of the signal and noise is no longer the main parameters for the design of this filter.
Property 2: The phase delay of this new filter is fixed. For the traditional digital filter, the phase delay changes with the change of input signal’s frequency. If the central point of the filtered data is taken as the result of the filter and the data length is $n$. The data in filter window is $\left[{x}_{1},{x}_{2},{x}_{3},\dots ,{x}_{n}\right]$, at next sampling time the data will change to $\left[{x}_{2},{x}_{3},{x}_{4},\dots ,{x}_{n+1}\right]$, when sampling point ${x}_{n+1}$ moves to the center of filter window, the data becomes $\left[{x}_{nn/21},{x}_{nn/2},{\dots ,{x}_{n},x}_{n+1},{x}_{n+2},\dots ,{x}_{n+n/2}\right]$. So the phase delay is $n/2\bullet t$, where $t$ is the sampling time. It’s clear that the phase delay is not related to the frequency of the input signal, and only related to data length $n$ and sampling time $t$.
Fig. 11. The decomposition result when the window length is 20
Fig. 12. The decomposition result when the window length is 50
4. The design of controller
4.1. The compensation of phase delay
The picture of practical object and the control block diagram of electrohydraulic servo system are shown in Fig. 13 and Fig. 14, respectively. The presence of the filter in feedback will have some influences on amplitude and phase. The phase delay of traditional filter like Butterworth filter and Chebyshev filter always grows larger while frequency becomes higher. But the phase delay of this new improved HHT filter is only related to the length of filtering window $n$, and maintains $n/\text{2}$ sampling periods delay. Based on this feature, the phase delay could be compensated in controller.
Fig. 13. The picture of the practical object
Fig. 14. The system block diagram
This electrohydraulic servo system is a position control system, and its transfer function can be obtained by using system identification. The form of transfer function is shown in Eq. (13):
in which $A\left({z}^{1}\right)=1+{a}_{1}{z}^{1}+a{z}^{2}+{a}_{3}{z}^{3}$, $B\left({z}^{1}\right)={b}_{0}+{b}_{1}{z}^{1}+{b}_{2}{z}^{2}+{b}_{3}{z}^{3}$. ${z}^{n/2}$ in Eq. (13) shows that there is a $n/2$ orders delay in system. To compensate this delay, the input $X\left(k\right)$ becomes $X\left(k+n/2\right)$, and then the change of output $Y\left(k\right)$ is shown in Eq. (14):
After compensation, the phase delay is only related to system itself.
4.2. Stability analysis of controller
The property 2 of this new filter means that a delay link is added to the system. The system’s open loop transfer function will be changed from $B\left({z}^{1}\right)/A\left({z}^{1}\right)$ to ${z}^{n/2}B\left({z}^{1}\right)/A\left({z}^{1}\right)$. If the input signal is $X\left(k+n/2\right)$, then the output signal $Y\left(k\right)$with noise $N\left(k\right)$is:
The noise after filtering becomes ${N}^{\text{'}}\left(k\right)$, ${N}^{\text{'}}\left(k\right)=N\left(k\right)e\left(k\right)$, $e\left(k\right)$ is the error between original noise and filtering result. Because the real filter cannot achieve the ideal effect, the error always exists. The system’s output after filtering is ${Y}^{\text{'}}\left(k\right)$, ${Y}^{\text{'}}\left(k\right)=\left(B\right({z}^{1})/A\left({z}^{1}\right))X\left(k\right)+e\left(k\right)$. So this new controller’s stability is only related to system itself, and does not introduce additional stability influence.
Fig. 15. The effect comparison chart at 2 Hz
Fig. 16. The effect comparison chart at 5 Hz
4.3. The control result
Apply this new controller in system and test the control performance at frequency 2 Hz, 5 Hz, 10 Hz, and 20 Hz. To illustrate the practicability and the superiority of this new controller, take Butterworth filter (Low pass filter, with pass band 020 Hz, pass band attenuation no more than 3 dB, the attenuation of 150 Hz at least 60 dB) as a comparison filter. The comparison results are shown in Figs. 1518, illustrating that the control error of this new filter is less than Butterworth filter especially at high frequency.
Fig. 17. The effect comparison chart at 10 Hz
Fig. 18. The effect comparison chart at 20 Hz
5. Conclusions
To design a new online filter based on HHT method, optimization of EMD method from four aspects are made and a better filtering performance is obtained. These optimizations are used to improve the effect of envelop curve which is the key of decomposition process. A new controller based on this new filter is also designed. Due to the selfadaptation characteristics of HHT method, the new controller’s phase delay is only related to the filtering window length $n$, and this delay could be compensated in controller. The control effect has been improved and more accurate control results are achieved.
References
 Wang Linhong, Wu Bo, Du Runsheng, et al. Nonlinear dynamic characteristics of moving hydraulic cylinder. Journal of Mechanical Engineering, Vol. 43, Issue 12, 2007, p. 1219, (in Chinese). [Search CrossRef]
 Huang N. E., Shen Z., Long S. R., et al. The empirical mode decomposition and the Hilbert spectrum for nonlinear and nonstationary time series analysis. Proceedings of the Royal Society of London Series A. Vol. 454, 1998, p. 903995. [Search CrossRef]
 Ramandeep Kaur, Vikramjit Singh Timefrequency domain characterization of stationary and nonstationary signals. International Journal for Research in Applied Science and Engineering Technology, Vol. 2, Issue 5, 2014, p. 438447. [Search CrossRef]
 Antonio H. C., Stephan H. Adaptive timefrequency analysis based of autoregressive modeling. Signal Processing, Vol. 91, Issue 4, 2011, p. 740749. [Search CrossRef]
 Fleureau J., Nunes J. C., et al. Turning tangent empirical mode decomposition: a framework for mono and multivariate signals. IEEE Trans on Signal Processing, Vol. 59, Issue 3, 2011, p. 13091316. [Search CrossRef]
 Marcus D., Torsten S. Performance and limitations of the HilbertHuang Transform (HHT) with an application to irregular water waves. Ocean Engineering, Vol. 31, 2004, p. 17831834. [Search CrossRef]
 Yang Peihong, Yue Lihai, Kang Lan, et al. Study on HHT based low frequency oscillation monitor for power system. Electrical Measurement and Instrumentation, Vol. 51, 21, p. 110114, (in Chinese). [Search CrossRef]
 Lu Zhimao, Jin Hui, Zhang Chunxiang, et al. Voice activity detection in complex environment based on HilbertHuang transform and order statistics filter. Journal of Electronics and Information Technology, Vol. 34, Issue 1, 2012, p. 213217, (in Chinese). [Search CrossRef]
 Park J. H., Lim H., Myung N. H. Analysis of jet engine modulation effect with extended HilbertHuang transform. Electronics Letters, Vol. 49, Issue 3, 2013, p. 215216. [Search CrossRef]
 Jie Huang, Jian Tang, Meijun Zhang, et al. An improved EMD based on cubic spline interpolation of extremum centers. Journal of Vibroengineering, Vol. 17, Issue 5, 2015, p. 23932409. [Search CrossRef]
 Li Zhao, Zhou Xiaojun, Xu Yun End effect treatment for EMD based on the period superposition extrapolation of mean generating function. Journal of Vibration and Shock, Vol. 32, Issue 15, 2013, p. 138143, (in Chinese). [Search CrossRef]
 Guo Mingwei, Ni Shihong, Zhu Jiahai, et al. HHT/EMD end extension method in vibration signal analysis. Journal of Vibration and Shock, Vol. 31, Issue 8, 2012, p. 6266, (in Chinese). [Search CrossRef]