Application of ESMD Method to Air-Sea Flux Investigation

The ESMD method can be seen as a new alternate of the well-known Hilbert-Huang transform (HHT) for non-steady data processing. It is good at finding the optimal adaptive global mean fitting curve, which is superior to the common least-square method and running-mean approach. Take the air-sea momentum flux investigation as an example, only when the non-turbulent wind components is well extracted, can the remainder signal be seen as actual oscillations caused by turbulence. With the aid of −5/3 power law for the turbulence, a mode-filtering approach based on ESMD decomposition is developed here. The test on observational data indicates that this approach is very feasible and it may greatly reduce the error caused by the non-turbulent components.


Introduction
The ESMD method is the abbreviation of "extreme-point symmetric mode decomposition" method [1] which can be seen as a new alternate of the well-known Hilbert-Huang transform (HHT) method [2], and it is anticipated as an effective usage in the fields of atmospheric and oceanic sciences, informatics, economics, ecology, medicine and seismology etc.Now the corresponding paper is available at http://arxiv.org/abs/1303.6540 with a protection of software by the national copyright administration of China [3,4].
Historically, in the field of data processing there are three typical methods: (1) The classical Fourier transform method, it is suitable for linear and stationary data; (2) The popular Wavelet transform method, it is suitable for linear and non-stationary data; (3) The hot HHT method based on the empirical mode decomposition (EMD), it is suitable for nonlinear and non-stationary data.At present, HHT is the most advanced method and is widely used.But there are also some shortcomings in it, such as: a) How to choose the stoppage criteria is a puzzling problem [5]; b) The residual component is too rough to reflect the whole evolutionary trend; c) The Hilbert-spectrum analyzing approach has some inherent defects.
Our ESMD method has not only remedied the shortcomings of HHT but also developed several new features as follows: (1) the outer envelop interpolation is changed to the inner extreme-point symmetric one; (2) the resi-dual component is optimized by the sifting times; (3) in view of the defect of all the integral transforms in timefrequency analysis, the traditional approach is abandoned and a data-based "direct interpolating" one is developed.
For a given data, its frequency analysis should be done on the oscillating part.Hence, to clear away the global mean curve is the first and foremost problem.In fact, only when the global mean curve is an optimal one, can the remainder signal be seen as actual oscillations caused by a series of wave fluctuations.In this sense, the ESMD method also offers a good adaptive approach for finding the optimal global mean fitting curve.It is superior to the commonly used least-square method (requires a priori function form) and running-mean approach (lacks of theoretical base and different choices of time-window and weight coefficients may result in different curves).This advantage is good for air-sea flux investigation.

On the Air-Sea Flux
The importance of air-sea interaction to the earth's climate is widely appreciated and this attracts the researchers to make more in situ flux observations at sea.Among all the air-sea flux measurement methods, the commonly used are the direct eddy-covariance method [6] and the inertial-dissipation method [7,8].The first one is based on the initial definition of the flux, which performs well on fixed platforms and it requires velocity corrections in case the platform is a moving one [9][10][11][12][13][14].The second one is based on the inertial-dissipation arguments of the tur-bulence and it is insensitive to the lower-frequency movements of the platform.From viewpoint of reliability the direct one is preferred, after all, the indirect one is associated with some uncertain parameters.
In terms of initial definition, not only the vector quantity (the momentum) but also the scalar quantity (the heat, moisture and CO 2 , etc.) relies on the vertical wind speed w.So in case w > 0 there should be an upward flux, and on the contrary there should be a downward flux.This is always true for the scalar quantities since the material concentrations are always positive.But it is not trivial for the momentum.Before calculating the flux, the signs of the horizontal wind velocities u and v should be dropped, since their signs only used to distinguish the horizontal directions.
Take the vertical transferring of the observed u-direction momentum as an example.Notice that the transport in the air-sea boundary layer is governed almost entirely by turbulence, the Reynolds decompositions ( ) ( ), ( ) ( ) u t u u t w t w w t are usually adopted.Here , u w and ( ), ( ) u t w t ′ ′ are the averages and fluctuations.In fact, only when the signal is a steady one can the fluctuations be seen as turbulent components.For the steady case, the total u-direction momentum flux is (2) Yet it is not always the case.The observed wind is usually an unsteady signal and the fluctuations may contain some non-turbulent components.So the adoption of the above formula may lead error to the air-sea flux.Efforts are done in [12][13][14] for it.In fact, the non-turbulent components ( ), ( ) u t w t can be seen as the synthesis of series lower-frequency modes extracted by ESMD decomposition.At this time, the total flux formula becomes

t w t u t u t w t w t u t w t u t w t u t w t u t w t
Here only the last term is due to the wind turbulence and it should be the real air-sea flux in statistics.Though the other three terms also contribute to the flux at the observational point, they are not recommended, after all, the non-turbulent flux is usually due to the horizontal unsteadiness.

Momentum Flux Investigation
The dataset used here (see Figure 1) came from an air-sea flux experiment on a moored ship in 2008, which is supported by the National High Technology Research and Development Program of China.Notice that the ship is almost fixed in the dock and the apparatus is installed almost vertically with 0, w ≈ the rotating correction on them is omitted here.Since the horizontal direction of the momentum is almost determined by the total means u and , v here we only care about its magnitude.To sim- plify the problem, if u and v do not change their signs at the chosen period, we can firstly take the absolute values of them before implementing the ESMD decomposition.
For the sign-sifting case as our data we can synthesize them first with The synthesized result is given in Figure 2.
Firstly, we take the Fourier spectrum of w (Figure 3).Notice that the energy spectrum of turbulence almost obeys −5/3 power law in the inertial sub-range, there are two referable approaches for filtering the non-turbulent components.The first one is a direct high-pass filtering approach used in [12,13], it relies on finding the critical frequency; the second one is a mode-filtering approach used in [14] which takes HHT as the basis, it relies on  determining the gap mode with the help of noise-assisted processing (called "ensemble EMD").Many flux calculating tests show that the first approach is very sensitive to the frequency truncation and a modification is needed on decomposing the wind eddies.In this sense, the second one is progressive.But the noise-assisted processing may destroy the signal and leads to a unbelievable result.Here we develop a direct mode-filtering approach based on ESMD.
By implementing the ESMD decomposition on w it yields the modes and frequency distributions in   To get rid of the reminder R and Modes 6 -8 it yields a Fourier spectrum in Figure 6.For this case the high frequency part accords well with the −5/3 power law.We note that the existence of lower frequency part here is very natural, after all, the Fourier spectrum is a result of linear transform.For this case, the synthesis of modes 1 -5 can be seen as turbulent part ( ). w t ′ Similarly, the turbulent part for U(t) can be also got by this approach.Take here the mean wind speed 4.978 m/ s U = at 8.8 m is transformed to that at 10 m according to the wind profile formula in [15] with an approximate drag coefficient 10 0.0024.

C =
The same magnitude of Equations ( 4) and ( 5) indicates the feasibility of our approach.By the way,  if the formula in forms of Equation ( 2) is directly used, it yields a flux 2 0.980 N/ m τ = which is too high to believe.

Figure 1 .
Figure 1.The observed wind components by HS-50 type of sonic anemometer with 20 Hz sampling rate at 8.8 m height.

Figure 2 .
Figure 2. The synthesized result of u and v.

Figure 4 and 5
separately.It follows from Figure3that the lower

Figure 3 .
Figure 3.The Fourier spectrum of w (the white line has a slope −5/3).

Figure 4 .
Figure 4.The ESMD decomposition of w with optimal 3 sifting times.
sign insures positive downward.Relatively, it follows from the common bulk formula that

Figure 5 .
Figure 5.The frequency distributions of the modes 4 -8 to w along the time.

Figure 6 .
Figure 6.The Fourier spectrum of modes 1 -5 to w (the