## 2 Materials and Methods
### MTV Normalization  
In the traditional TV method, the sensitivity to smoothing parameters is significantly high due to nonlinearity caused by data fitting, regularization terms, and convergence. As a result, the LSRTM problem becomes highly unstable, and its convergence cannot be guaranteed. The objective function of the LSRTM-TV is as follows:  
\[S^{\epsilon}(\mathbf{m})=\underset{\mathbf{m}}{\text{min}}\Big{\{}\|d_{obs}- L\mathbf{m}\|_{2}^{2}+\lambda\|\mathbf{m}\|_{TV}\Big{\}} \tag{10}\]  
where \(\|\cdot\|_{TV}\) is the TV regularization operator as follows:  
\[\|\mathbf{m}\|_{TV}=TV(\mathbf{m})=\int_{\Omega}|\
abla\mathbf{m}|d\Omega \tag{11}\]  
where \(\Omega\) is the imaging region, \(|\cdot|\) is the absolute value sign, and \(\
abla\) is the gradient operator. Since the derivative of the TV operator is non-continuous, it is guaranteed to be differentiable by the following approximation:  
\[TV_{\delta}(\mathbf{m})=\int_{\Omega}\sqrt{|\
abla\mathbf{m}|^{2}+\delta^{2}}d\Omega \tag{12}\]  
Hence, we introduce a new LSRTM method with MTV regularization to improve the imaging resolution and mitigate artifacts. The objective function of the LSRTM-MTV method is given as follows:  
\[S^{\epsilon}(\mathbf{m},\mathbf{u})=\underset{\mathbf{m},\mathbf{u}}{\text{ min}}\Big{\{}\|d_{obs}-L\mathbf{m}\|_{2}^{2}+\lambda_{1}\|\mathbf{m}-\mathbf{u}\|_{2}^ {2}+\lambda_{2}\|\mathbf{u}\|_{TV}\Big{\}} \tag{13}\]  
where, \(\lambda_{1}\) and \(\lambda_{2}\) are both positive regularization parameters. Compared with a conventional \(TV\) regularization term, Equation (13) contains a new variable \(\mathbf{u}\) and an additional regularization term. The objective function can be written as:  
\[S^{\epsilon}(\mathbf{m},\mathbf{u})=\underset{\mathbf{u}}{\text{min}}\Big{\{} \underset{\mathbf{m}}{\text{min}}\Big{\{}\|d_{obs}-L\mathbf{m}\|_{2}^{2}+ \lambda_{1}\|\mathbf{m}-\mathbf{u}\|_{2}^{2}\Big{\}}+\lambda_{2}\|\mathbf{u}\| _{TV}\Big{\}} \tag{14}\]As can be seen from Equation (14), the regularization parameter \(\lambda_{1}\) controls the trade-off between the data misfit term and the Tikhonov regularization term, and \(\lambda_{2}\) is employed here to balance the amount of interface preservation in L8 RTM. The alternating-minimization algorithm is used to solve the dual minimization problem [31]. By using an initializer \(\mathbf{u}^{(0)}=\mathbf{m}^{(0)}\) to solve Equation (14), the solutions of two minimization sub-problems can be obtained:  
\[\mathbf{m}^{(k)}=\underset{\mathbf{m}}{\text{argmin}}\{S_{1}^{\varepsilon}( \mathbf{m})\}=\underset{\mathbf{m}}{\text{argmin}}\Big{\{}\|d_{\text{obs}}-L \mathbf{m}\|_{2}^{2}+\lambda_{1}\|\mathbf{m}-\mathbf{u}^{(k-1)}\|_{2}^{2}\Big{\}} \tag{15}\]  
\[\mathbf{u}^{(k)}=\underset{\mathbf{m}}{\text{argmin}}\{S_{2}^{\varepsilon}( \mathbf{m})\}=\underset{\mathbf{u}}{\text{argmin}}\Big{\{}\|\mathbf{m}^{(k)}- \mathbf{u}\|_{2}^{2}+\lambda_{2}\|\mathbf{u}\|_{TV}\Big{\}} \tag{16}\]  
where \(k\) is the iteration step. Noting that the two sub-problems correspond to different parts, the first is to solve for \(\mathbf{m}^{(k)}\) by using L8 RTM with the Tikhonov regularization [32] and prior image \(\mathbf{u}^{(k-1)}\); while the second is to solve for \(\mathbf{u}^{(k)}\) using the L2-TV minimization to preserve the outline sharpness of interfaces within the L8 RTM image \(\mathbf{m}^{(k)}\), as well, it also denoises the image. In this paper, the Limited-memory Broyden-Fletcher-Goldfarb-Shanno (L-BFGS) algorithm [33] is employed to solve for m(k) in the sub-problems (15), while the Split-Bregman method [34] is utilized to solve for u(k) in the sub-problems (16). The iterative formula for \(\mathbf{m}\) is as follows:  
\[\mathbf{m}^{(k+1)}=\mathbf{m}^{(k)}+\mathbf{a}^{(k)}\mathbf{p}^{(k)} \tag{17}\]  
where \(\mathbf{p}^{(k)}\) is the updating direction and \(\mathbf{a}^{(k)}\) is the iteration step size. Initially, the L-BFGS method is employed to obtain parameter \(\mathbf{p}^{(k)}\) using the calculated function \(\mathbf{g}\) based on Equation (8).