Impulse noise treatment in magnetotelluric inversion

Hugo Hidalgo-Silva; Enrique Gómez-Treviño

doi:10.1515/geo-2020-0225

Article Open Access

Impulse noise treatment in magnetotelluric inversion

Hugo Hidalgo-Silva and Enrique Gómez-Treviño

Published/Copyright: February 3, 2021

Published by

Become an author with De Gruyter Brill

Submit Manuscript Author Information

From the journal Open Geosciences Volume 13 Issue 1

Abstract

The problem of model recovering in the presence of impulse noise on the data is considered for the magnetotelluric (MT) inverse problem. The application of total variation regularization along with L1-norm penalized data fitting (TVL1) is the usual approach for the impulse noise treatment in image recovery. This combination works poorly when a high level of impulse noise is present on the data. A nonconvex operator named smoothly clipped absolute deviation (TVSCAD) was recently applied to the image recovery problem. This operator is solved using a sequence of TVL1 equivalent problems, providing a significant improvement over TVL1. In practice, TVSCAD requires the selection of several parameters, a task that can be very difficult to attain. A more simple approach to the presence of impulse noise in data is presented here. A nonconvex function is also considered in the data fitness operator, along with the total variation regularization operator. The nonconvex operator is solved by following a half-quadratic procedure of minimization. Results are presented for synthetic and also for field data, assessing the proposed algorithm’s capacity in model recovering under the influence of impulse noise on data for the MT problem.

Keywords: impulse noise; magnetotelluric data; geosounding

1 Introduction

Electromagnetic sounding methods have been applied for a long time to investigate the interior of the Earth from depths of a few meters to hundreds of kilometers. The magnetotelluric (MT) sounding method is one of the most popular today, and it is aimed at estimating the electrical conductivity distribution at depth associated with subsurface geological structures. It has been successfully applied in shallow exploration problems and to infer the Earth’s deep electrical conductivity distribution, temperature regimes, and geological structure. The method uses measurements of natural electric and magnetic fields recorded at the surface of the Earth. If it is assumed that the underground structure is one-dimensional (1D), the inverse problem consists of recovering the vertical conductivity distribution from the ratio of surface electric and magnetic field measurements taken at different frequencies. One of the main characteristics of the MT inverse problem is its nonlinearity. Several methods have been developed to cope with this problem [1]: linearization, asymptotic behavior, and exact methods. Constable et al. [2] consider smooth models using a procedure called Occam’s inversion. The method consists of Tikhonov’s regularization applied to the linearized equations of the 1D MT problem. A fundamental characteristic of Occam’s models is that the resulting subsurface conductivity distribution is a continuous function. This is very convenient from a mathematical point of view [3]. Due to the roughness penalizer inclusion of Constable et al. [2], the developed model tends to smear discontinuities, an undesirable feature in most situations. In an effort to generate models with sharp boundaries, several techniques have been applied to this problem: piecewise continuous formulations as in ref. [4,5], l ₁ norm minimization [6,7,8], total variation penalizers in ref. [9] and [10], among others. An important requirement for the regularizer is to allow the recovery of edges and smooth the homogeneous parts. As is well known, total variation [11] is now the standard approach to meet this requirement. TV consists on applying a convex regularizer, ensuring existence and uniqueness of a solution. In the case of regular functions, the total variation of u is the L 1 norm of the gradient, i.e., TV( u ) = ∫ Ω | ∇ u | d x .

There are two basic methods for MT transfer function estimation from field data: the single-site method and the remote reference method. In the first method, one or two biased estimates at each frequency are obtained for each measuring direction [12]. In the second method, only one unbiased estimate can be obtained at one site simultaneously using a magnetic field measurement at a remote site [13,14]. The latter can be used effectively only for uncorrelated noise between local and remote sites. The remote reference method is a standard procedure in MT data acquisition because is the only way to remove bias. The single-site method is still used when a second instrument is not available or is faulty. The MT linear relationship of the impedance tensor is given by

(1) E = Z H + N

where E is the local electric field, Z is the MT transfer function tensor, H is the magnetic field, and N is the noise. The impedance Z is either estimated by the single source or remote reference method. The noise N accounts for higher dimensional sources and short duration non-stationary and instrumental noise. N vanishes for the strict zero wavenumber or plane wave model. Once Z is estimated, the residuals can be computed by substituting their values in equation (1). Residuals are given by the difference between the observed E field and the predicted E field, obtained by substituting the estimated Z in relation (1). The use of classical spectral analysis together with least-squares regression is warranted if data follow a stationary and Gaussian model. In this case, it is also common to assume that residuals follow a multivariate normal probability distribution. If this is satisfied, the maximum likelihood estimates are obtained. In practice, most data display gross departures, or outliers, from the simple regression model [15]. The main causes are geomagnetic phenomena, thunderstorms, anthropogenic contribution, and instrumental problems. Outliers usually appear as a fraction of the useful observations, with distinct characteristics from the rest of the sample. The phenomenon is independent of the structure of the bulk of geomagnetic field observations. In this paper, we do not consider the robust estimation of the MT transfer functions, the treatment of MT data contaminated with large variations is the main purpose here.

Results are presented for the MT problem, but can also be applied to other types of electromagnetic soundings.

2 MT inversion

Tikhonov’s regularization method is the standard technique applied to obtain models of the subsurface conductivity distribution from electric or electromagnetic measurements. The model proposed for the conductivity distribution (m) is obtained by minimizing the functional

(2) U T ( m ) = μ ∥ F ( m ) − d ∥ ℋ 2 + P ( m ) ,

with F : ℳ ⊃ D ( F ) → ℋ representing the direct functional, applied to an element on the model space ℳ , a Banach space and returning a member of the Hilbert ℋ space of data. P ( m ) = ∥ ∇ m ∥ ℳ 2 is for Tikhonov’s maximum smoothness functional the usual approach. Due to the roughness penalizer inclusion, the model developed by Tikhonov’s algorithm tends to smear discontinuities, a feature that may be undesirable in situations where different rock types simply face each other without a fuzzy transition.

In an effort to generate models with sharp boundaries, several techniques have been applied to this problem: piecewise continuous formulations as in ref. [4,5], l 1 norm minimization [6,7,8], total variation penalizers in ref. [9,10], among others. An important requirement for the regularizer is to allow the recovery of edges and smooth the homogeneous parts. As is well known, Total Variation [11] is now the standard approach to meet this requirement. Considering P ( m ) = TV ( m ) in (2), a more realistic model is recovered, but the resulting model can be severely affected if the Gaussian’s assumption of data noise is not accomplished. The algorithm solving (2) with P ( m ) = TV ( m ) and l 2 norm for fitting will be called TVL2. The model developed by solving TVL2 is observed in Figure 1a for the synthetic data presented in Figure 1b, obtained by solving the forward equations for a four-layered model.

$Figure 1 Developed model and fitness to apparent resistivity when solving μ ∥ d − F ( m ) ∥ 2 + TV ( m ) \mu \parallel d-F(m){\parallel }^{2}+\text{TV}(m) with Gaussian noise: (a) developed model and (b) fitness to ρ a.$

Figure 1

Developed model and fitness to apparent resistivity when solving μ ∥ d − F ( m ) ∥ 2 + TV ( m ) with Gaussian noise: (a) developed model and (b) fitness to ρ _a.

In Figure 2a, a model is obtained by the same algorithm solving (2) but inserting a large variation on the first amplitude data. Notice the effect on the resulting model.

$Figure 2 Recovered model and fitness to apparent resistivity for μ ∥ d − F ( m ) ∥ 2 + TV ( m ) \mu \parallel d-F(m){\parallel }^{2}+\text{TV}(m) considering Gaussian noise with impulses at first two and last data: (a) model and (b) magnitude fitness.$

Figure 2

Recovered model and fitness to apparent resistivity for μ ∥ d − F ( m ) ∥ 2 + TV ( m ) considering Gaussian noise with impulses at first two and last data: (a) model and (b) magnitude fitness.

For the case of impulse noise, a preferred model is TVL1, replacing the l 2 norm in the data fitness term with l 1 , and considering TV as regularizer [16]:

(3) U TVL 1 ( m ) = μ ∥ F ( m ) − d ∥ 1 + ∑ i ∥ D i m ∥ ,

where D i m denotes a local first-order finite-difference of the ith element of m. It has been observed that l 1 -norm as a fidelity term removes the influence of abnormal values in the case of image denoising [17]. However, the TVL1 model may be severely affected in the case of high noise levels. The reason is that all data, both corrupted and noise-free, are equally penalized in the data fitting term [18], leading to great difficulty for balancing regularization and data fitting. A model obtained by TVL1 for the same data as before is observed in Figure 3a.

$Figure 3 Recovered model and fitness to apparent resistivity for μ ∥ d − F ( m ) ∥ 1 + TV ( m ) \mu \parallel d-F(m){\parallel }_{1}+\text{TV}(m) considering Gaussian noise with impulses at first two and last data: (a) model and (b) magnitude fitness.$

Figure 3

Recovered model and fitness to apparent resistivity for μ ∥ d − F ( m ) ∥ 1 + TV ( m ) considering Gaussian noise with impulses at first two and last data: (a) model and (b) magnitude fitness.

Nikolova [19] noticed that the solutions of the TVL1 model substantially deviate from both the data-acquisition and the prior models. Gu et al. [18] proposed to combine TV regularization with nonconvex smoothly clipped absolute deviation penalty for data fitting called TVSCAD. A difference of convex functions is adopted to solve the non-convex TVSCAD.

Cai et al. [20] proposed a two-phase approach to restore images corrupted by blur and impulse noise called CTVL1. The first phase is to detect possibly corrupted pixels, then in the second phase, they restore images using only the noncorrupted pixels. The adaptive correction procedure obtains an initial estimator of the image by TVL1 and then apply a corrective step by blurring the initially estimated image [20]. The connections among TVSCAD and CTVL1 are observed by Gu et al. [18]. They notice that the main difference is in the function considered in the correction step of CTVL1. They claim to obtain comparable results to TVL1. An important issue of TVSCAD is the selection of SCAD function parameters. As observed by Gu et al., they are problem-dependent, but the best choices are very difficult to determine.

Recently, Zhang et al. [21] proposed a nonconvex approach to image restoration contaminated with impulse noise. They consider a nonconvex data fitting term and TV for regularization operator. The nonconvex functions used are an exponential type and the Geman function P ( t ) = t t + β . They proposed a proximal linearized minimization algorithm for the nonconvex part and the alternate direction method of multipliers (ADMM) for the total variation part. They sow convergence for the case of a concave function P ( m ) .

In this paper, we also consider applying a nonconvex function to the impulse noise treatment. The nonconvex part is solved using a half-quadratic minimization algorithm proposed by Geman and Yang [22]. The total variation regularized is also considered with the ADMM method. The algorithm is easy to implement and the convergence properties were established by Charbonnier et al. for a nonconvex function satisfying the edge-preserving properties [23]. The application is to the MT nonlinear inverse problem. As mentioned earlier, other algorithms (TVSCAD, CTVL1) present some difficulties in parameter selection, and others like [21] require a concave function to converge.

2.1 Nonconvex function applied to noise treatment

A smooth, nonconvex regularization algorithm was proposed by Charbonnier et al. [23], considering the following operator:

(4) P EPR = ∑ j ϕ [ ( D x m ) j ] + ∑ j ϕ [ ( D z m ) j ] ,

where ( D x m ) j is the finite-difference implementation of first-order derivative respect to x at point j, ( D z m ) j is the corresponding derivative respect to z , and ϕ is the potential function. This operator has been considered on the vision community under the name of “edge-preserving regularization” (EPR) [23] and applied to the resistivity inverse problem by Hidalgo et al. [24]. Charbonnier et al. developed some properties that a potential function has to comply in order to obtain an EPR operator; some functions satisfying the properties are presented by them. A particularly important nonconvex potential function initially proposed by Geman and McClure [25] for image restoration is

(5) ϕ ( x ) = x 2 x 2 + β 2 ,

with β > 0. This function does not satisfy assumptions required by Zhang et al. for convergence [21] because it is not concave as x → 0 .

Charbonnier et al. apply an iterative algorithm developed by Geman and Yang [22] called “half-quadratic regularization” for the minimization of a modified version of the original functional (2) with P = P _EPR. This algorithm is basically applying a convex conjugate to replace the original nonconvex function with a new one defined over an extended domain. Charbonnier et al. showed that for a potential function satisfying several conditions, the sequence generated by this algorithm is convergent.

In this paper, we propose to apply the operator to the fitness term in (2), solving

(6) U TVNC ( m ) = 1 2 N μ M ∑ i ϕ ( r M i ) + μ P ∑ i ϕ ( r P i ) + P TV ( m ) ,

with ϕ(⋅) applied to the magnitude r M i = ( d M i − log ρ a ( m , ω i ) ) , and phase r P i = ( d P i − p ( m ⋅ ω i ) ) residuals of fitness to measurements at frequency ω i . ρ a is the apparent resistivity at the surface ρ a ( m , ω i ) = | Z 1 | 2 ω μ o , with Z 1 the impedance of the superior layer on a layered media as described on ref. [4], and μ o is the magnetic permeability of free space. p ( m . ω i ) = arctan Im ( Z 1 ) Re ( Z 1 ) is the phase of the impedance at frequency ω _i. A first-order Taylor approximation for F(m) is considered for the minimization, F ( m + Δ m ) = F ( m ) + J ( m ) Δ m + ϵ , with J ( m ) = ∂ F ∂ m , assuming a small step Δm and neglecting ϵ . Minimization is realized by applying Charbonnier et al. iterative algorithm “half-quadratic regularization” to the fitness term, and considering a Bregman split algorithm to the regularization operator.

The algorithm is aimed to minimize the functional

(7) U TVNC 2 ( m , b M , b P ) = 1 2 N μ M ∑ i b M i r M i 2 + ψ ( b M i ) + μ P ∑ i b P i r P i 2 + ψ ( b P i ) + ∑ j | ν j | + γ 2 ∑ j ( m j − m j − 1 − ν j − d j ′ ) 2

with two auxiliary variables b M i , b P i for magnitude ( r M ) and phase ( r P ) residuals, respectively. Also, two regularization parameters μ _M and μ _P are included to consider the different residual magnitudes. d′ is an auxiliary variable used to implement the Bregman iteration, updated after solving (7). ψ is the convex dual function satisfying ϕ ( t ) = inf w ( w t 2 + ψ ( w ) ) . The minimization procedure is based on alternating minimizations over m , b , and ν . Minimization with respect to m is quadratic, using coordinate descent, the solution for rth element is

(8) m r = 1 N ∑ j d ˆ j − ∑ l ≠ r w j l m l ) ( w j r μ ˆ j + γ n 1 ( m r ) 1 N ∑ j w j r 2 μ ˆ r + 2 γ

with w i j = ∂ m i ∂ F j , d ˆ i = d i − F i ( m ) + ∑ j w i j m j k , μ ˆ i = μ M b M i or μ ˆ i = μ P b P i , depending on the corresponding data, and n 1 ( m r ) = m r − 1 + ν r + d r + m r + 1 + ν r + 1 + d r + 1 is the first-order neighborhood of m r .

The solution with respect to v is given by a shrinkage operator,

(9) v r k + 1 = shrink ( D ( m r k + 1 ) − ( d ′ r ) k , γ )

with shrink ( u , β ) = max ( | u | − β , 0 ) u | u | , and minimum on b are given by

(10) b M j = ϕ ′ ( r M j ) 2 r M j

(11) b P j = ϕ ′ ( r P j ) 2 r P j .

With ϕ′ the derivative of ϕ.

The Geman–McClure function is applied here, but several other functions can also be considered, as observed in ref. [23].

Parameter selection is not as complicated as in other methods; μ M and μ P can be selected by some of the regularization parameter selection methods available [26]. The penalty parameter γ is obtained by trial and error, and as noticed by Getreuer [27], the algorithm is not so sensitive to the fair value of γ. Parameter β is selected as a scaling term according to the model’s dynamic range.

The algorithm is presented in Table 1, minimization is stopped when the model does not change in some min value.

Table 1

Iterative algorithm for (7)

Algorithm

Input: d

Initialization: k = 0 , m 0 = 0 , d 0 = 0

repeat

evaluate F ( m k ) , J ( m k ) , W = F ( m k ) − m k J ( m k )

d k + 1 ← d + ( d k − W m k )

m k + 1 ← arg min m U TVNC 2 ( m , b M k , b P k , v k )

ν k + 1 ← arg min v U TVNC 2 ( m k + 1 , b M k , b P k , v )

update b M k + 1 , b P k + 1 with (10) and (11)

k = k + 1

until ( ∥ m k − m k − 1 ∥ < min )

Output m k

It should be noticed that, as the model evolves from a homogeneous single layer to a multilayered system, the residuals will be large in the first steps. Because of that the potential function should be adapted, to consider those large residuals at the first steps. To implement that parameter β is initialized to a value so that a quadratic error fitness function is approximated for a few (ten) steps. Then, the desired value of β is applied to obtain a non-convex fitness function. The algorithm parameters are observed in Table 2.

Table 2

Algorithm parameters

Algorithm	μ M	μ P	γ	RMS
TVL2	10	10	10	8.4
TVL1	100	10	10	8.22
TVNC2	100	10	10	8.63

TVL1 and TVL2 were stopped when ∥ m k − m k − 1 ∥ was less than 0.04 and 0.01, respectively. All algorithms were asked to stop when min was 0.01, but TVL1 would not converge for the parameters shown in Table 2. TVL1 is implemented using ADMM for both, TV and L1 minimization, it uses other parameter γ for the L1 part, its value was also 10.

For TVL2, both regularization parameters were the same because the standard deviation is used as a weight for the measurements; a standard deviation of 0.05 for magnitude and 0.5 degrees for phase were considered in the example. The standard deviation affects the selection of μ, which is the reason for a smaller μ in that algorithm. The developed model by TVNC2 and the corresponding fitness are shown in Figure 4(a and b). It can be observed that the recovery of the model improves significantly over the previous cases. Notice that the resistivity of the third layer is underestimated. This is not so much a drawback of the algorithm but an intrinsic feature of the MT method, which detects much better good conductors than good resistors. This is because in good resistors, the induced currents are smaller than in good conductors. Consequently, the former does not alter significantly the electric and magnetic fields on the surface.

$Figure 4 Recovered model and fitness to apparent resistivity for μ ∑ i ϕ ( r i ) + TV ( m ) \mu {\sum }_{i}\phi ({r}_{i})+\text{TV}(m) after applying algorithm TVNC2, considering Gaussian noise with impulses at first two and last data: (a) model and (b) magnitude fitness.$

Figure 4

Recovered model and fitness to apparent resistivity for μ ∑ i ϕ ( r i ) + TV ( m ) after applying algorithm TVNC2, considering Gaussian noise with impulses at first two and last data: (a) model and (b) magnitude fitness.

2.2 Application to field data

Figure 5 presents a model developed by TVNC2 for field data of site PC5008 in the COPROD2 data set [28]. The data are the “uncorrected” original data provided by Jones [28]. The field data are observed in Figure 6(a and b). A large deviation can be observed in the first two lectures.

Figure 5

Model obtained by algorithm TVNC2 for site PC5008 of uncorrected data of COPROD2.

Figure 6

Original, uncorrected data fitness for site PC5004 of COPROD2 by TVNC2 algorithm: (a) log apparent resistivity fitness and (b) phase fitness.

It should be noticed that TVL2 does not converge for this uncorrected data. A correction was made by Jones [28], obtaining the data observed in Figure 7(b and c). The model developed by TVL2 is observed in Figure 7a. Notice the similarity with the one obtained by the TVNC2 algorithm when using the uncorrected data.

Figure 7

Model obtained by TVL2 and data fitness for the corrected data of site PC5008 of the COPROD2 data set: (a) model, (b) magnitude fitness, and (c) phase fitness.

3 Conclusion

The MT sounding data are usually contaminated with large errors. The robust treatment of MT data has been generally applied in the transfer function components determination. In this paper, the case of large errors on the resulting apparent resistivity and phase data is considered.

A comparison is made among the usual least-squares inversion, a robust technique based on the L1 norm applied on the fitness term, and a nonconvex operator proposed in the image restoration some time ago for image noise filtering.

The application of a nonconvex operator for the fitness term and TV for regularization for impulse noise treatment on MT data is proposed. The inverse problem is solved by using a half-quadratic procedure for the nonconvex operator, along with a split Bregman procedure, ensuring convergence. The proposed algorithm parameters are easy to select.

Results are presented for data obtained from a synthetic model, observing the capability of the proposed method in model recovering when large deviations are inserted on the data. Application to field COPROD2 data is presented too, assessing the capability of the algorithm to recover a model of raw, uncorrected data.

Acknowledgments

The authors would like to thank P. J. Savage of PanCanadian Petroleum Limited of Calgary, Alberta, for making the raw (uncorrected) data available to the Geological Survey of Canada.

Author contributions: HHS developed and implemented the algorithms, and EGT assisted with the analysis and interpretation of results.

References

[1] Whittall KP, Oldenburg DW. Inversion of magnetotelluric data for a one-dimensional conductivity. Tulsa: Society of Exploration Geophysicists; 1992 Jan 1.10.1190/1.9781560802419Search in Google Scholar

[2] Constable SC, Parker RL, Constable CG. Occam’s inversion: A practical algorithm for generating smooth models from electromagnetic sounding data. Geophysics. 1987 Mar;52(3):289–300.10.1190/1.1442303Search in Google Scholar

[3] Gómez-Treviño E. Nonlinear integral equations for electromagnetic inverse problems. Geophysics. 1987 Sep;52(9):1297–302.10.1190/1.1442390Search in Google Scholar

[4] Hidalgo H, Marroquin JL, Gómez-Treviño E. Piecewise smooth models for electromagnetic inverse problems. IEEE Trans Geosci Remote Sens. 1998 Mar;36(2):556–61.10.1109/36.662738Search in Google Scholar

[5] Varentsov IM. A general approach to the magnetotelluric data inversion in a piecewise-continuous medium. Izv. Phys Solid Earth C/C Fiz Zemli-Rossiiskaia Acad. Nauk. 2002 Nov 1;38(11):913–34.Search in Google Scholar

[6] Esparza FJ, Gómez-Treviño E. Electromagnetic sounding in the resistive limit and the Backus-Gilbert method for estimating averages. Geoexploration. 1987 Dec 1;24(6):441–54.10.1016/0016-7142(87)90013-5Search in Google Scholar

[7] Esparza FJ, Gómez-Treviño E. 1-D inversion of resistivity and induced polarization data for the least number of layers. Geophysics. 1997 Nov 1;62(6):1724–9.10.1190/1.1444272Search in Google Scholar

[8] Hidalgo H, Gómez-Treviño E, Pérez-Flores MA. Linear programs for the reconstruction of 2d images from geophysical electromagnetic measurements. Subsurface Sens Technol Appl. 2004 Apr 1;5(2):79–96.10.1023/B:SSTA.0000034544.25382.8fSearch in Google Scholar

[9] Dobson DC. Recovery of blocky images in electrical impedance tomography. In Inverse problems in medical imaging and nondestructive testing. Vienna: Springer; 1997. p. 43–64.10.1007/978-3-7091-6521-8_5Search in Google Scholar

[10] Vogel CR, Oman ME. Fast, robust total variation-based reconstruction of noisy, blurred images. IEEE Trans Image Process. 1998 Jun;7(6):813–24.10.1109/83.679423Search in Google Scholar

[11] Rudin LI, Osher S, Fatemi E. Nonlinear total variation based noise removal algorithms. Phys D Nonlinear Phenom. 1992 Nov 1;60(1–4):259–68.10.1016/0167-2789(92)90242-FSearch in Google Scholar

[12] Sims WE, Bostick Jr, FX, Smith HW. The estimation of magnetotelluric impedance tensor elements from measured data. Geophysics. 1971 Oct;36(5):938–42.10.1190/1.1440225Search in Google Scholar

[13] Gamble TD, Goubau WM, Clarke J. Magnetotellurics with a remote magnetic reference. Geophysics. 1979 Jan;44(1):53–68.10.1190/1.1440923Search in Google Scholar

[14] Clarke J, Gamble TD, Goubau WM, Koch RH, Miracky R. Remote-reference magnetotellurics: equipment and procedures. Geophys Prospecting. 1983 Jan 1;31(1):149–70.10.1111/j.1365-2478.1983.tb01047.xSearch in Google Scholar

[15] Trad DO, Travassos JM. Wavelet filtering of magnetotelluric data. Geophysics. 2000 Mar;65(2):482–91.10.1190/1.1444742Search in Google Scholar

[16] Bai M, Zhang X, Shao Q. Adaptive correction procedure for TVL1 image deblurring under impulse noise. Inverse Probl. 2016 Jun 16;32(8):085004.10.1088/0266-5611/32/8/085004Search in Google Scholar

[17] Nikolova M. A variational approach to remove outliers and impulse noise. J Math Imaging Vis. 2004 Jan 1;20(1–2):99–120.10.1023/B:JMIV.0000011920.58935.9cSearch in Google Scholar

[18] Gu G, Jiang S, Yang J. A TVSCAD approach for image deblurring with impulsive noise. Inverse Probl. 2017 Nov 10;33(12):125008.10.1088/1361-6420/aa9383Search in Google Scholar

[19] Nikolova M. Model distortions in Bayesian MAP reconstruction. Inverse Probl & Imaging. 2007;1(2):399.10.3934/ipi.2007.1.399Search in Google Scholar

[20] Cai JF, Chan RH, Nikolova M. Two-phase approach for deblurring images corrupted by impulse plus Gaussian noise. Inverse Probl Imaging. 2008;2(2):187.10.3934/ipi.2008.2.187Search in Google Scholar

[21] Zhang X, Bai M, Ng MK. Nonconvex-TV based image restoration with impulse noise removal. SIAM J Imaging Sci. 2017;10(3):1627–67.10.1137/16M1076034Search in Google Scholar

[22] Geman D, Yang C. Nonlinear image recovery with half-quadratic regularization. IEEE Trans Image Process. 1995 Jul;4(7):932–46.10.1109/83.392335Search in Google Scholar PubMed

[23] Charbonnier P, Blanc-Féraud L, Aubert G, Barlaud M. Deterministic edge-preserving regularization in computed imaging. IEEE Trans Image Process. 1997 Feb;6(2):298–311.10.1109/83.551699Search in Google Scholar PubMed

[24] Hidalgo H, Gómez-Treviño E, Marroquin JL, Esparza FJ. Piecewise continuous models for resistivity soundings. IEEE Trans Geosci Remote Sens. 2001 Dec;39(12):2725–8.10.1109/36.975007Search in Google Scholar

[25] Geman S, McClure D. Statistical methods for tomographic image reconstruction. Bull Int Stat Inst. 1987;52:5–21. In Proceedings of the 46th Session of the International Statistical Institute.Search in Google Scholar

[26] Vogel CR. Computational methods for inverse problems. Philadelphia: Society for Industrial and Applied Mathematics; 2002 Jan 1.10.1137/1.9780898717570Search in Google Scholar

[27] Getreuer P. Rudin–Osher–Fatemi total variation denoising using split Bregman. Image Process Line. 2012 May 19;2:74–95.10.5201/ipol.2012.g-tvdSearch in Google Scholar

[28] Jones AG. The COPROD2 dataset: Tectonic setting, recorded MT data, and comparison of models. J Geomag Geoelec. 1993 Sep 20;45(9):933–55.10.5636/jgg.45.933Search in Google Scholar

Received: 2020-03-23

Revised: 2021-01-15

Accepted: 2021-01-16

Published Online: 2021-02-03

This work is licensed under the Creative Commons Attribution 4.0 International License.

Articles in the same Issue

https://doi.org/10.1515/geo-2020-0225

Keywords for this article

impulse noise; magnetotelluric data; geosounding

Creative Commons

BY 4.0