## 2 Simulation Methods of Non-Stationary Stochastic Processes
### Spectral Representation Using the GHW Method  
The concept of HW and GHW were first proposed by [PERSON] [27] in 1994. Compared with other wavelets, HW has a complete set of orthogonal basis functions. However, its resolution is not enough for the analysis of some fine structures. Then, GHW with a higher resolution was proposed in due course. As defined by [PERSON], in the time-domain, the GHW function has the following structure:  
\[w_{(m,n),k}\left(t-\frac{k}{n-m}\right)=\frac{\left\{\exp\left(i\cdot 2n\pi \left(t-\frac{k}{n-m}\right)\right)-\exp\left(i\cdot 2m\pi\left(t-\frac{k}{n-m} \right)\right)\right\}}{i\cdot 2\pi(n-m)\left(t-\frac{k}{n-m}\right)} \tag{7}\]  
in which \(m,n\) are different real numbers that denote different wavelet levels, \((n-m)\) denotes the bandwidth and \(n>m,k\) is the translation transformation factor and \(k/(n-m)\) represents the time step, and \(i^{2}=-1\) is the imaginary unit. In order to intuitively see the influence of different frequency bandwidths and time shift factors on GHW, the schematic diagram of equally spaced GHW with different frequency bandwidths is given below. Let \(\Delta=n-m\) denote the equally spaced bandwidth. Figures 1-3 show the real and imaginary parts of GHW for the case that \(n=7\), \(m=6\), \(\Delta=1\), \(k=0\); \(n=8\), \(m=6\), \(\Delta=2\), \(k=0\); and \(n=7\), \(m=6\), \(\Delta=1\), \(k=1\), respectively. Figure 1: Real and imaginary parts of GHW at \(n=7\), \(m=6\), \(\Delta=1\), \(k=0\). Figure 3: Real and imaginary parts of GHW at \(n=7\), \(m=6\), \(\Delta=1\), \(k=1\). In fact, from the expression of Equation (7), the real part is an even function and the imaginary part is an odd function. Comparing Figures 1 and 2, it can be concluded that, as the bandwidth increases, the GHW waveform becomes narrower. The energy of GHW is mainly focused on the time interval \([(k/(n-m))-1,(k/(n-m))+1]\) and the value range of real and imaginary parts is the same no matter how wide the frequency bands are. Comparing Figures 1 and 3, it can be found that the position of the main energy is related to the time-varying translation factor \(k\). There is no obvious difference between the real part and the imaginary part at the same scale \(m\), \(n\), \(k\), except that the real part reaches the maximum value and the imaginary part is 0 when \(t=k/(n-m)\). The corresponding Fourier transform of Equation (7) gives the expression of GHW transform in the frequency domain as follows:  
\[W_{(m,n),k}(\omega)=\left\{ \begin{array}{cc}\{1/(n-m)\cdot 2\pi\}\exp\{-i \omega k/(n-m)\}&2\pi m\leq\omega\ \leq 2\pi n,\\ 0&\text{elsewhere}.\end{array} \right. \tag{8}\]  
in which \(2\pi n\), \(2\pi n\) denote different angular frequencies. Equation (8) shows that \(W_{(m,n),k}\) (\(\omega\)) is equal to zero outside the interval \([2\pi n\), \(2\pi n]\) and is constant within this interval. The GHW coefficients \(c_{(m,n),k}\) and \(\widetilde{c}_{(m,n),k}\) are also known as GHW transform and can be defined as integral of the product of a time process and GHW function. \[ \begin{array}{c}c_{(m,n),k}=(n-m)\int_{-\infty}^{+\infty}f(t)\cdot w_{(m,n), k}^{*}\Big{(}t-\frac{k}{n-m}\Big{)}\text{d}t\\ \widetilde{c}_{(m,n),k}=(n-m)\int_{-\infty}^{+\infty}f(t)\cdot w_{(m,n),k} \Big{(}t-\frac{k}{n-m}\Big{)}\text{d}t\end{array}  \tag{9}\]  
As derived by [PERSON] and [PERSON] [29], a stochastic process \(f\) (\(t\)) can be represented by the GHW function \(w_{(m,n),k}\) and coefficient \(c_{(m,n),k}\) as follows:  
\[f^{\text{GHW}}(t)=\sum_{m,n}\sum_{k}\sqrt{2\pi\cdot(n-m)\cdot S_{(m,n),k}} \cdot w_{(m,n),k}\cdot c_{(m,n),k} \tag{10}\]  
in which \(S_{(m,n),k}\) represents the EPSD at scale \(m\),\(n\) and translation \(k\).