**3. Methods**

The expressions (11) and (12) contain Sommerfeld integrals, which is an infinite integral with the highly oscillatory and slow-decaying kernel. In this paper, the partial closed form of Sommerfeld integral is derived, and ESPRT is applied to extract DCIM, while DE rules are used for the computation of the finite integration.

#### *3.1. Partial Closed Form Expression*

The Sommerfeld integration in Equations (11) and (12) can be written in the following form:

$$I = \int\_0^\infty \frac{m}{m\_i} \mathbf{g}(m\_i) e^{-m\_i|z - z\_s|} f\_0(mr) dm = \int\_0^\infty \Gamma\_m f\_0(mr) dm \tag{17}$$

The integration is then divided into two parts,

$$\int\_{-0}^{\infty} \Gamma\_{m} l\_{0}(mr)dm = \int\_{0}^{p} \Gamma\_{m} l\_{0}(mr)dm + \int\_{-p}^{\infty} \Gamma\_{m} l\_{0}(mr)dm\tag{18}$$

where *p* is a reasonably selected integral breakpoint. The second integral on the right side of the Equation (18) can be asymptotically approximated and then be written as

$$\int\_{-p}^{\infty} \Gamma\_m f\_0(mr) dm \approx \int\_{-0}^{\infty} \Gamma\_m^{\infty} f\_0(mr) dm - \int\_{-0}^{p} \Gamma\_m^{\infty} f\_0(mr) dm \tag{19}$$

On substituting (19) in (18), we ge<sup>t</sup>

$$I \approx \int\_{-0}^{p} (\Gamma\_m - \Gamma\_m^{\infty}) f\_0(mr) dm + \int\_{-0}^{\infty} \Gamma\_m^{\infty} f\_0(mr) dm \tag{20}$$

The kernel function *g*(*mi*) in (17) can be approximated by an exponential function,

$$\Gamma\_{m}^{\infty} = g(m\_i) \cdot \frac{m}{m\_i} e^{-m\_i|z - z'|} = \sum\_{l=1}^{pb} a(l) \exp\left[b(l)m\_i\right] \cdot \frac{m}{m\_i} e^{-m\_i|z - z'|} \tag{21}$$

where *pb* is the number of exponentials used for approximation. In this paper, the coefficients *a* (*l*), *b* (*l*) are solved by the ESPRIT algorithm, which will be discussed in the next section.

From the Sommerfeld identity expressed by formula (5), the closed form of the second integral can be obtained,

$$\int\_{-0}^{\infty} \Gamma\_m^{\infty} l\_0(mr) dm = \sum\_{l=1}^{pb} a \ (l) \frac{\exp(-jk\_l R\_l)}{R\_l} \tag{22}$$

where *Rl* = 3*r*2 + [*b* (*l*) − |*z* − *zs*|]2, *r* = 3(*x* − *xs*)<sup>2</sup> + (*y* − *ys*)2, (*xs*, *ys*, *zs*) is the location of the source.

The first part of the *g*(*mi*) function usually contains singularity, and the tail is smooth and decays fast. Therefore, if the appropriate *p* value is selected, the approximate fitting will be very accurate by avoiding the singular value in the front part. The finite integral with singularity in (20) can be evaluated directly by the numerical integration method. We choose DE quadrature rules here for integration with the advantage of dealing with singular points and high precision.
