Hilbert transform using a robust geostatistical method

In this paper, we introduced an efficient inversion method for Hilbert transform calculation which can be able to eliminate the outlier noise. The Most Frequent Value method (MFV) developed by Steiner merged with an inversion-based Fourier transform to introduce a powerful Fourier transform. The Fourier transform process (IRLS-FT) ability to noise overthrow efficiency and refusal to outliers make it an applicable method in the field of seismic data processing. In the first part of the study, we introduced the Hilbert transform stand on a efficient inversion, after that as an example we obtain the absolute value of the analytical signal which can be used as an attribute gauge. The method depends on a dual inversion, first we obtain the Fourier spectrum of the time signal via inversion, after that, the spectrum calculated via transformation of Hilbert transforms into time range using a efficient inversion. Steiner Weights is used later and calculated using the Iterative Reweighting Least Squares (IRLS) method (efficient inverse Fourier transform). Hermite functions in a series expansion are used to discretize the spectrum of the signal in time. These expansion coefficients are the unknowns in this case. The test procedure was made on a Ricker wavelet signal loaded with Cauchy distribution noise to test the new Hilbert transform. The method shows very good resistance to outlier noises better than the conventional (DFT) method.


Introduction
Seismic attributes play an important role in seismic data processing and interpretation. They can be physical attributes (further classified as pre-stack or post-stack ones having sub-classes of instantenous or wavelet attributes) or geometrical attributes used to obtain certain information about dip, azimuth and discontinuity (further classifications are also exist). They are used for detecting fractures or faults and geological lineaments which are below the seismic resolution for providing important structural and stratigraphic information in three dimensions.
Amplitude, phase, and frequency are the main characteristics elements of seismic traces. Extracting instantaneous and local values of these parameters is a useful way of revealing hidden features. In addition to the basic attributes, weighted averaging of instantaneous seismic attributes is beneficial. This technique is most useful for reducing spikes and noise variations [1].
On another important side, it's become an important tool to geoscientists for reservoir characterization and quality control. Since the paper published by Taner et al. [2], the subject area has expanded and significantly developed. Today, we can talk about a wide range of attributes, divide them into physical and geometric attributes, classify procedures according to whether they are applied before or after processing, interpreted on one or more seismic channels, etc. In developing and IOP Conf. Series: Earth and Environmental Science 942 (2021) 012029 IOP Publishing doi:10.1088/1755-1315/942/1/012029 2 applying measurement and data processing procedures, it is important to manage and improve the signal-to-noise relationship as much as possible.
Fourier transform is often a priority in the creation of attribute stations. Szegedi and Dobroka [3] proposed a robust Fourier transform process (IRLS-FT) built on an inversion basis. In this paper, we demonstrate that the method works effectively in suppressing outlier noise and can repair the signalto-noise relationship by up to an order of magnitude. As a first step, we present the production of the Hilbert transform as part of a transformed robust/resistant inversion, which plays a decisive role in the definition of the complex channel. The reflection strength attribute was introduced in the article which represents the total energy strength of the trace and is independent of phase due to which it is more suitable to analyze the amplitude anomalies in detail.

The analytical signal
The starting point for calculating the basic attribute stations is the creation of an analytical signal (analytical or complex channel). The concept of analytical signal in data processing was introduced by the Nobel Prize laureate Hungarian physicist Dénes Gábor [4]. His ambition was to use the powerful mathematical tools of quantum mechanics in signal processing (using square-integrable complex functions as elements of the so-called Hilbert space). To do this he introduced the analytical signal where ( ) is the time signal and is the Hilbert -transform of the time signal. The analytical signal or complex trace was first discussed in details by Bracewell [5] and Taner [2] defining the most important seismic attributes: the reflection strength, the instantaneous phase and frequency, etc. The concept and its use is well-discussed is seismic literature. According to Eq. (2), the Hilbert transform is generated as a convolution of the time signal u(t) with function t 1 π − . In the frequency domain, this relation can be written as where F denotes the Fourier transform. As giving the Hilbert transform as It can be seen that when calculating the Hilbert transform the Fourier transform and its inverse are to be used. It is well known that the traditional DFT and IDFT algorithms are noise-sensitive in the case of non-Gaussian noise in the data set. To define a robust procedure an inversion-based Fourier transform method (IRLS-FT) was introduced by Dobroka et al. [6]. In order to use the noise-rejection capacity of the method, we can use IRLS-FT in calculating ( ) (6). A full inversion-based method can be produced when the

Reflection strength attribute and its noise sensitivity
Once the analytical signal (6) is known, the attributes can be produced. The reflection strength (instantaneous amplitude, paving) is further examined as an example, which is the absolute value of the analytical signal The noise sensitivity of this attribute is illustrated using the Ricker wavelet [7] shown in Figure 1 (in blue). The time series shows a 10 Hz wave packet localized at 0.1 sec with a sampling interval of 0.005 sec in the [-1,1] (sec) domain. The data set (I) was generated by adding Gaussian noise (with a standard deviation of 0.0025 and zero mean) to the noise-free data. A data system (II) with outlier errors was created by generating noise following the Cauchy distribution of scale parameter ε=0.04. Figures 1-3 show the reflection strength channel calculated by the conventional process for noise-free inputs laden with Gauss noise (I) and Cauchy noise (II). Hilbert's transformed functions were carried out using the discrete Fourier transform process (DFT) and its inverse (IDFT).  It can be seen that the reflection strength gauge calculated based on the data system containing outliers II is particularly noisy. To characterize noise sensitivity, we introduce distance in the data space, which is d(I) = 0.0240 for data set I and for the data set II the distance d(II) = 0.0444 (N being the number of the samples). Figure 3 justifies the development of a Hilbert transform process more resistant to outlier noise. Since Hilbert's transformed training is based on the traditional Fourier transform (DFT), it is obvious that the inversion-based Fourier transform procedure (IRLS-FT) utilizing Steiner weights [3] is used to solve the task.

The inversion-based Fourier transform (IRLS-FT)
The Fourier transform links the time domain of the signal registration and the frequency domain of the signal test according to the following formulas: The frequency spectrum U(ω) is a Fourier transform of the time signal u(t), which is usually a continuous function of the complex value. When using series expansion based discretization, the spectrum is written in terms of a suitably chosen system of basis functions: where denotes the complex expansion coefficients, Ψ (ω) denotes the n-th known basis function, and M is the total number of the basis functions. If the Fourier transform is understood as an overdetermined inverse problem, we must first designate the direct problem, which is the inverse Fourier transform, defined in the case of the k-th measurement datum as Combining the above equations, the calculated (theoretical) data can be determined according to the following system of linear equations: is the Jacobi matrix. The deviation vector of measured and calculated data is given by The expression contains the expansion coefficients, which are determined by minimizing some norm of the deviation vector. Knowing the expansion coefficients, the spectrum can be determined at any frequency as

The robust generation of Hilbert transform
To produce the Hilbert transform we need to know the U(ω) spectrum of the signal. To improve the signal-to-noise ratio, the Fourier transform is performed by the IRLS-FT method, and then the spectrum is multiplied by the function -j sgn(ω). We then return to the time domain by an inverse Fourier transform. This can be done by applying IDFT and also by using an inversion-based inverse Fourier transform. The latter procedure can be performed by robust inversion defined here. The starting point is the expression of the inverse Fourier transform The direct problem is given by the formula of the Fourier transform, where the time function u(t) is discretized in the form of a series expansion After substitution, the following formula is obtained for the k-th sampling element of the spectrum is the Jacobi matrix, the elements of which can be thought of as Fourier transforms of the basis function system. The calculation of the complex integral in the formula can be avoided by choosing the basic functions of series expansion from the eigenfunctions of the Fourier transform.

Choosing basis function system -the Hermite functions
The Hermite functions have been chosen as basis functions because in the case of inverse problems it is advised to use complete, orthonormal function systems to reduce the number of unknown parameters and to improve the stability of the inversion procedure. The Hermite functions require proper scaling because in geophysical applications the frequency covers a wide range. The scaled Hermite polynomials can be calculated using the Rodriguez formula where α is a scaling factor. The first and second Hermite polynomials are The scaled Hermite polynomials fulfil the condition of orthogonality In the formula, δnm denotes the Kronecker-delta symbol. The scaled Hermite function can be defined as where hn(ω,α) denotes the n-th scaled Hermite polynomial. Thus the introduced Hermite functions are orthonormal To get a fast and simple formula for the calculation of Jacobi's matrix a special feature of the Hermite functions was used. As we can see earlier, Gk,n is the inverse Fourier transform of the basis function. This is the main reason for selecting the Hermite functions because their non-scaled versions are eigenfunctions of the Fourier transform [9] and can be written for the inverse Fourier transform Jacobi's matrix can be expressed by the scaled Hermite functions Hn(ω,α) as Introducing the following notations Jacobi's matrix can be written again based on Eq. (15) This is a very important result because Jacobi's matrix can be derived quickly without integration. The discretized form of the spectrum can be written according to Eq. (3) using the scaled Hermite function system, where the expansion coefficients Bn are defined (including the expression of the Jacobi's matrix) in the frame of the over-determined inverse problem. This inverse problem is highly overdetermined because the number of measurement data is much higher than that of the parameters (N>M). In the case of the Least Squares method (LSQ), the L2 norm of the deviation vector is minimized The well-known normal equation can be derived from the above condition in the following form After this, we can estimate the complex series expansion coefficients (35) With this, the real and imaginary part of the Fourier spectrum can be calculated at any frequency by Eq. (3) . The LSQ method gives optimal results in the case of a Gaussian distributed data set.
To make the Fourier transform more robust an Iteratively Reweighted Least Squares (IRLS [5]) method using Cauchy weights was implemented. To achieve the optimal values of the unknown parameters (Bn) the following weighted norm is minimized where the weighs are defined as and the k-th element of the deviation vector is The scale parameter ε of the Cauchy distribution is not a priori given because the data residuals change from iteration to iteration [10]. The weighted norm gives reliable results for inverse problems even if the measured data set contains outliers. The j-th iteration step of the normal equation is This iteration is repeated until a proper stop criterion is met. Finally, the Fourier spectrum can be calculated at any frequency by using Eq. (10).

Numerical tests
In Figures 1, 2, and 3, examples were shown of noise-free and the absolute value of the analytical signal (reflection strength) calculated with a Ricker wavelet laden with Gaussian and Cauchy distribution noise. In Figures 4 and 5, the absolute value of the analytical signal (the reflection strength, as the first seismic attribute) generated by the robust inversion-based Fourier transforms of the same input signals is illustrated. Of course, our method for noise-free input provides the same result as the Fourier transform using the traditional DFT process. In the case of input data system I (Gauss noise), the inversion-based Hilbert transform in Figure 4 shows that it is less noisy than the data in Figure 2 produced by the traditional method. The slight improvement is characterized by a distance of d(I)=0.0183 in the data space. However, significant improvement is reflected in the processing of data system II (Cauchy noise) using the robust inversion method defined with Steiner weights as shown in Figure 5. Here, compared to Figure 3, we can see an almost complete suppression of the effect of the outlier data, which is caused by the distance in the data space d(II) =0.0033. There is nearly an order of magnitude difference between the data space distances obtained by the traditional (d(II) =0.0444) and inversion-based (d(II) =0.0033) Fourier transform processes.   Figure 5. Result after applying the new algorithm on the noisy signal using Cauchy distribution noise (data set II), we can see that the result is very effective on those type of noises.
The above results confirm a significant improvement in signal/noise when using a robust inversionbased Fourier transform and an inverse Fourier transform built on a robust inversion basis to produce the Hilbert transform.

Conclusions
The new IRLS inversion-based Fourier transform method has been demonstrated to achieve high noise resistance to excessively noisy data. In this paper, we used this procedure in the calculation of the Hilbert transform. The algorithm presented is based on a dual inversion, on the one hand, we apply IRLS-FT to the Fourier transform and, on the other hand, after the production of the transformed spectrum of Hilbert, the inverse Fourier transform is also calculated using an inversion-based robust/resistant process. Although dual inversions require a significant calculation time compared to traditional DFT/IDFT transform, the numerical example presented shows that we can achieve a high degree of improvement with the inversion-based Hilbert transform process. Given the current capacity and speed of computers, it is likely that in some practical cases the extra computing time will be tolerable.